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ABSTRACT 


We report the discovery of the unusually bright long-duration gamma-ray burst (GRB), 
GRB 221009A, as observed by the Neil Gehrels Swift Observatory (Swift), Monitor of All-sky X- 
ray Image (MAXI), and Neutron Star Interior Composition Explorer Mission (NICER). This energetic 
GRB was located relatively nearby (z — 0.151), allowing for sustained observations of the afterglow. 
The large X-ray luminosity and low Galactic latitude (b = 4.3?) make GRB 221009A a powerful probe 
of dust in the Milky Way. Using echo tomography we map the line-of-sight dust distribution and 
find evidence for significant column densities at large distances (2 10kpc). We present analysis of 
the light curves and spectra at X-ray and UV/optical wavelengths, and find that the X-ray afterglow 
of GRB 2210094 is more than an order of magnitude brighter at Ty + 4.5 ks than any previous GRB 
observed by Swift. In its rest frame GRB 221009A is at the high end of the afterglow luminosity 
distribution, but not uniquely so. In a simulation of randomly generated bursts, only 1 in 10^ long 
GRBs were as energetic as GRB 2210094; such a large Ey iso implies a narrow jet structure, but the 
afterglow light curve is inconsistent with simple top-hat jet models. Using the sample of Swift GRBs 
with redshifts, we estimate that GRBs as energetic and nearby as GRB 221009A occur at a rate of < 
1 per 1000 yr — making this a truly remarkable opportunity unlikely to be repeated in our lifetime. 


Keywords: Gamma-ray bursts (629)— High energy astrophysics (739)— Interstellar dust (836) 


1. INTRODUCTION 


Massive stars exhibit a broad continuum of proper- 
ties in their terminal explosions. At one end are the 
cosmological long-duration gamma-ray bursts (GRBs), 
capable of coupling tremendous energies (> 10°! erg) 
to highly collimated ejecta with a bulk Lorentz factor 
Io Z 100 (Piran 2004). On the other end, the energy 
budget of most stripped-envelope core-collapse super- 
novae is dominated by the (quasi)-isotropic supernova 
emission, with photospheric velocities of tens of thou- 
sands of kms ! (e.g., Liu et al. 2016). Intermediate 
between these two extremes lies the growing class of 
low-luminosity GRBs and relativistic supernovae (e.g., 
Margutti et al. 2014). Typified by the prototypical 
low-luminosity GRB 980425 associated with SN 1998bw 
(Galama et al. 1998; Kulkarni et al. 1998), these sources 
couple several orders of magnitude less energy to their 
moderately relativistic ejecta (~ 1075 erg), and lack the 
high degree of collimation of cosmological GRBs. 

In the nearby universe (z S 0.3), where high-energy 
facilities are sensitive to low-luminosity GRBs, high- 
luminosity GRBs are exceedingly rare due to their much 
lower volumetric rate. Yet an energetic GRB observed 
within this volume could produce unprecedented bright- 
ness. Here we report the discovery of GRB 2210094, 
an extremely luminous GRB (Kann & Agui Fernandez 
2022) in our cosmic backyard. 

On 2022 October 9 at 14:10:17 UT, the Burst Alert 
Telescope (BAT; Barthelmy et al. 2005), onboard the 
Neil Gehrels Swift Observatory (Swift; Gehrels et al. 
2004), triggered twice in rapid succession on a new cos- 
mic source in the constellation Sagitta. Following its 
automated burst response, Swift promptly slewed to 


the location of the first trigger, detecting a bright tran- 
sient seen with both the Swift X-Ray Telescope (XRT; 
Burrows et al. 2005) and the Ultraviolet/Optical Tele- 
scope (UVOT; Roming et al. 2005). Due to the rarity of 
such a repeated BAT image trigger and the proximity of 
the source to the Galactic plane (b — 4.3?), it was ini- 
tially classified as a new Galactic X-ray and optical tran- 
sient, and therefore was designated Swift J1931.14-1946 
(Dichiara et al. 2022). Monitor of All-sky X-ray Image 
(MAXI; Matsuoka et al. 2009) reported the detection of 
bright X-ray emission from this location shortly there- 
after (Negoro et al. 2022). 

Subsequently the Gamma-ray Burst Monitor (GBM; 
Meegan et al. 2009) onboard Fermi reported the de- 
tection of an exceptionally bright long-duration GRB 
7:55 min prior to the initial BAT trigger, with a consis- 
tent localization (Veres et al. 2022). Due to issues in 
receiving data, the automated classification and local- 
ization notices associated with this onboard GBM trig- 
ger were not distributed to the world. However, the 
Fermi team rapidly communicated the existence of the 
GBM trigger, and its spatial coincidence with the dou- 
ble BAT trigger, to the Swift team. This spatial associ- 
ation, together with analysis of prompt XRT data that 
showed a smooth GRB-like power-law decline, led to the 
conclusion that Swift J1931.14-1946 was in fact a GRB, 
GRB 221009A (Kennea et al. 2022). For the first time 
in the ~18 years since launch, BAT triggered not on the 
GRB prompt emission, but instead on the bright high- 
energy afterglow when GRB 221009A entered the field 
of view. 

The unusual brightness of GRB 221009A prompted 
widespread follow-up at multiple wavelengths. Addi- 


tional X-ray detections were reported by the Neutron 
Star Interior Composition Explorer (NICER; Iwakiri 
et al. 2022) and the Nuclear Spectroscopic Telescope 
Array (NuSTAR; Brethauer et al. 2022). The Large 
High Altitude Air Shower Observatory reported detect- 
ing photons up to 18 TeV (Huang et al. 2022). Changes 
in the strength of signals propagated by radio transmit- 
ters were recorded at the time of the GBM trigger as 
the photons from the GRB ionized Earth's atmosphere 
(Schnoor et al. 2022; Guha & Nicholson 2022). Spectro- 
scopic observations of the afterglow and the host galaxy 
provided a redshift of z — 0.151 (de Ugarte Postigo et al. 
2022; Castro-Tirado et al. 2022; Izzo et al. 2022), corre- 
sponding to a distance of 749.3 Mpc. 

'This paper is organized as follows: $2 contains analy- 
sis of the observations taken by Swift (BAT, XRT, and 
UVOT), MAXI, and NICER; in 83, we present anal- 
ysis of the dust scattering echo, broadband spectrum, 
and light curve of the burst afterglow; we discuss how 
GRB 221009A compares to other GRBs in 84, investi- 
gate the astrophysical rate of similar events, and the 
nature of energetic GRBs; in $5, we present our conclu- 
sions. We show that due to the combination of proxim- 
ity and large (but not unprecedented) intrinsic luminos- 
ity GRB 221009A has a much brighter X-ray afterglow 
than previously observed Swift GRBs, and such lumi- 
nous nearby events are extremely rare occurrences. 

For this paper, we assume a cosmology with Hp = 
67.36, Qm = 0.3153, Q, = 0.6847 (Planck Collaboration 
et al. 2020), and z — 0.151 unless otherwise stated. We 
adopt the GBM trigger time as the burst onset, i.e., To 
= 13:16:59.99 UTC on 2022 October 9. Magnitudes are 
reported on the Vega system, and uncertainties are given 
at a 9076 confidence interval (unless otherwise noted). 


2. OBSERVATIONS 
2.1. Swift Burst Alert Telescope 


Figure 1 shows the BAT raw light curves (i.e., 
not background-subtracted), summed over all detec- 
tors, from To — 500s to To + 5000s. The location of 
GRB 221009A was occulted by the Earth until Tp + 
1870s. At ~To + 1100s, the overall count rate began 
to rise due to increased particle background as Swift 
approached the South Atlantic Anomaly (SAA). From 
To+1317-2183s, BAT data collection was disabled as 
Swift transited the SAA. The count rate remained ele- 
vated while exiting the SAA until the spacecraft slew 
beginning at To + 2496s, when we attribute the lin- 
early decaying enhanced count rate to emission from 
GRB 221009A. However, the source location was not in 
the coded portion of the BAT field of view until a slew 
beginning at To + 3095s. 


3 


Finally after this slew completed, BAT triggered on 
GRB 221009A, at 14:10:18 (To + 3199s; trigger ID 
1126853) and 14:17:06 (To + 3607s; trigger ID 1126854) 
UTC. The event data from these two triggers cover a 
time range from To + 2960s to To + 4570s. The mask- 
weighted light curve shows steadily declining emission 
present when the burst location came into the BAT 
field of view at Tọ + 3173s, and extending beyond the 
available event data range. The time-averaged spectrum 
from To + 3302s to To + 4538s is best fit by a simple 
power-law model, with T = 2.08 + 0.03. The fluence in 
the 15-150 keV band is (7.4 X: 0.1) x 1075 ergcm ?. 

'The smooth temporal evolution observed by BAT, to- 
gether with the large temporal offset between the BAT 
and GBM triggers, indicate that BAT triggered on the 
afterglow of GRB 2210094. This marks the first such oc- 
currence of a BAT afterglow trigger in the 18 years of 
Swift operations. 

Given the exceptionally bright afterglow, we searched 
for even later emission in the BAT survey mode data us- 
ing the BatAnalysis python package (Parsotan et al., in 
prep.). In individual pointings of survey mode data the 
afterglow was detected until 2022 October 9 21:55:38 UT 
(To + 31ks). We attempted to fit the spectra of each 
survey dataset with a power-law (CFLUX*PO) model in 
XSPEC (Arnaud 1996) to obtain fluxes and photon in- 
dices. Non-detections were then analyzed to obtain 5c 
upper limits following the procedure outlined in Laha 
et al. (2022). 

At later times, the survey data were binned daily and 
then mosaiced together. The spectrum from each mo- 
saiced image was fitted, and the flux and photon indices 
derived similar to the procedure above. The results of 
this analysis are plotted in the top panel of Figure 2, 
and listed in Appendix A. 


2.2. Swift X-Ray Telescope 


The XRT began observing GRB221009A at 
14:13:09 UT (To + 3370s, 170s after the first BAT 
trigger) and located a bright afterglow in the initial 
0.1-s image-mode exposure, after which observations 
began in Windowed Timing (WT) mode. The initial 
WT count rate was 910 40cts Í (all XRT count rates 
are corrected for the effects of pile-up and hot columns, 
see e.g., Evans et al. 2007); making it more than an 
order of magnitude brighter at this time — in observed 
flux - than any other GRB observed by XRT. Due to 
the high count rate, the XRT remained initially in WT 
mode, in which only one-dimensional spatial informa- 
tion is collected. Significant structures are present in 
the 1-D spatial profiles, with a clear excess compared to 
the expected point spread function (PSF), which were 
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Figure 1. BAT raw light curve from the rate data. Panel (a) shows the raw light curve in the 15-350 keV band with a time 
binning of 1.6s. This light curve was made from BAT quad-rate data, which records continuous count rates from all active 
detectors in four different energy ranges (15-25, 25-50, 50-100, and 100-350 keV). The time period when the GRB was occulted 
by the Earth (within 69? of the Earth center) is marked in green. Spacecraft slew times are marked in gray. The two BAT 
trigger times are marked by red lines. Panel (b) shows the partial coding fraction as a function of time. A partial coding fraction 
of zero indicates that the GRB was outside of the BAT coded field of view, and a value of one indicates that the source is in 


the highest sensitivity region of the BAT coded field of view. 


evolving with time (see Appendix B). This resembles 
the behavior expected when dust clouds in our Galaxy 
scatter X-rays from the GRB prompt emission, that 
were not initially traveling towards Earth, back into our 
line of sight. 

At To 4- 89 ks, the GRB had faded sufficiently that the 
XRT automatically switched to Photon Counting (PC) 
mode. The 2-D image from this observation, shown 
in Figure 3, confirmed the presence of a complex se- 
ries of expanding, bright rings associated with a dust- 
scattering echo (see also Tiengo et al. 2022), the prop- 
erties of which are discussed in § 3.1. 

'The presence of scattered emission complicates the 
data analysis, as it can contribute events to the regions 
over which source and background counts are accumu- 
lated. Both the intensity and spectrum of the rings 
are spatially variable, thus the selected ‘background’ 
region may in fact not be representative of the back- 
ground within the source region. Ás a result, the au- 
tomated XRT analysis (Evans et al. 2007, 2009) needs 
some modification. For PC mode this is relatively sim- 
ple, as we have full 2-D imaging and the innermost rings 
were reasonably well separated from the GRB itself. We 
restricted the source region to a radius of 20 pixels (47"), 
and manually defined à background region as free from 


dust contamination as possible (identified from a stacked 
image of all PC mode data). 

For WT mode this is more problematic since the CCD 
is read out in columns, rather than pixels. After in- 
vestigating the variation in the echo contribution to 
W'[-mode data, we modified the default WT extrac- 
tion regions to minimize its impact (see Appendix B). 
We found that the WT flux is accurate to ~6%. We fur- 
ther verified this by extracting a light curve using only 
data above 4keV, where dust scattering becomes less 
efficient, given the roughly v-? dependence of the scat- 
tering cross section; only minimal changes in the light 
curve shape are observed. On the other hand, the dust 
has more significant impact on the WT spectra (Ap- 
pendix B), resulting in an increase in the best-fit pho- 
ton index up to ~7% and a significant increase to the 
best-fit absorption column (up to ~27%). 

We modified the settings for the automated light 
curve! to limit the source extraction region to 20 pix- 
els radius, and to use the PC-mode background re- 
gion (described in the Appendix B); due to the unusual 
brightness of the source we also increased the number 


1 https:/ /www.swift.ac.uk/xrt curves/01126853/. 
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Figure 2. Combined X-ray and UV/optical light curve of GRB 2210094. The upper panel shows observed X-ray flux from 
MAXI, NICER, XRT (all 0.3-10.0 keV) and BAT (14-195 keV). BAT data are taken from detections of the afterglow in Survey 
data; note that the final BAT data point combines all observations integrated over 1 day. The dotted line shows a broken 
power-law fit to the MAXI, NICER, and XRT data. The inset figure shows the first MAXI, BAT and XRT data, with the 
black dashed line indicating a fit to the XRT data alone, and red dashed line a fit between the first MAXI point and the first 
XRT detection. The lower panel shows 7 filter optical and UV data from UVOT as obtained after the subtraction of late-time 
template images. The dashed line shows a power-law fit to the white band data. 


of counts per bin by a factor of ten compared to the 
default (see Evans et al. 2007 for details). We fit the 
time-averaged PC mode spectrum using XSPEC with a 
power-law model and two absorption components. The 
first of these was a TBABS fixed to the Galactic value of 
5.38x10?! em ? (Willingale et al. 2013); the second was 
a ZTBABS with Ng free and redshift fixed at 0.151. The 
fit yields Ng = (1.4+ 0.4) x 102 cm~?, with a photon 
index [ = 1.843 0.2. We used these results to convert the 
light curve into observed 0.3-10 keV flux. The resulting 
flux measurements are shown in Figure 2. 


To investigate possible spectral evolution, we ex- 
tracted a series of time-resolved spectra using the ‘Add 
time-sliced spectrum’ option on the UKSSDC website? 
(Evans et al. 2009), which ensures that the modifications 
to the default processing parameters, described above, 
are employed. From the W'T-mode data we created one 
spectrum from the first spacecraft orbit (T94-3.4—4.5 ks), 
and then spectra in 10ks chunks from 704-10—50 ks. 


? https:/ /www.swift.ac.uk/xrt.spectra/01126853/. 
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Figure 3. XRT PC-mode 0.8 — 5.0 keV images from ObsIDs 01126853004 (top left), 01126853006 (top right), 01126853008/09 
(bottom left) and 01126853010/11 (bottom right), illustrating the echo expansion with time. The echo event radial positions 
were scaled by t^^? to each observation mid-time (given in the upper left corner of each image in days since GBM trigger) to 


counter the halo expansion within each observation. 


For PC mode, we extracted spectra over the intervals 
To+67—200 ks, 200—400 ks, 400-1000 ks and > 1000 ks. 
For WT mode we used only grade 0 events, as is rec- 
ommended for very absorbed objects?. We fitted these 
spectra independently with the model defined above; the 
results are shown in Table 1. The photon index ranges 
between ~ 1.6 and ~ 1.9 with the harder profiles ob- 
served during the first block of observations (from 3.4 


3 https: //www.swift.ac.uk/analysis/xrt /digest_cal.php#abs 


to 4.5ks) and around 1-2 days (from 68 to 175ks) after 
Ty. 


2.3. Swift Ultraviolet/Optical Telescope 


'The UVOT began settled observations of the field of 
GRB 2210094 at 2022 October 9 14:13:17 UT, To+3.4 ks 
(179s after the first BAT trigger; Dichiara et al. 2022; 
Kuin et al. 2022). A nearby star, 5" away from the GRB, 
contaminates the photometry when using the standard 
circular aperture of 5" radius. In order to minimize the 


Table 1. Spectral parameters determined from fitting a power-law model to the time- 
resolved X-ray data from Swift, MAXI and NICER. 


T — To Exposure Intrinsic Ny Photon Fit stat. / Observatory / 
(ks) (ks) (102em ?) index d.o.f Detector 
2.5 0.047 — 1.75+0.09 — 307/322 | MAXI/GSC 
3.4 — 4.5 1.2 1.28 +0.04 1.6140.02 1249 / 959 XRT/WT 
8.0 0.047 — LOTTO e 145/166 — MAXI/GSC 
10.0 — 20.0 0.5 1.390.09 1.8940.05 859/840 — XRT/WT 
14.0 0.029 — 1.890.16 — 26.9/35 ^| NICER/XTI 
19.5 0.062 — 1.844 0.15  23.6/36 ^ NICER/XTI 
20.0 — 30.0 2.9 1.3640.05 1.84490; 995/922 — XRT/WT 
30.0 — 40.0 Ly 1.253-0.07 1.85+0.04 800/840 . XRT/WT 
40.0 — 50.0 2.8 1.20+0.06 — 1.8440.03 875/873 . XRT/WT 
52.6 — 98.0 2.1 — 1.85+0.07 131/153 NICER/XTI 
68 — 175 17 lr 1.65+0.08 — 655 / 727 XRT/PC 
103.0 — 198.0 6.5 — 1.792 +0.06 144/168 NICER/XTI 
203.3 — 360.5 9.8 — Ll 94/94  NICER/XTI 
216 — 394 12 1.264013 1.80 +0.10 632 / 665 XRT/PC 
404 — 954 26 1.20+0.18  1.92+0.10 501 / 622 XRT/PC 
505.0 — 961.6 2.7 — 1.82+0.19 — 58.1/67 ^ NICER/XTI 
1005 — 6285 0.24 1.004015 — 1.91 0.09 723 / 668 XRT/PC 


NoTE—NICER/XTI and MAXI/GSC data are fit in the energy range 4— 10 keV and 4— 
20 keV, respectively, to avoid contamination from unresolved dust echo contamination, 


with intrinsic Ny fixed to 1.29x107? em *. 


XRT spectra are fit in the range 0.3 — 


10 keV, allowing for measurement of intrinsic Ny. Note that fit statistics are C-stat 
(Cash 1979) for MAXI and Swift data, NICER are x”. Galactic Nu is fixed to a value 


of 5.38x10?! cm? (Willingale et al. 2013). 


star's contribution to the measurement of the afterglow, 
we use a 2.5" aperture to extract source counts. To be 
consistent with the UVOT calibration, these count rates 
were then corrected to 5" using the curve of growth con- 
tained in the calibration files. Background counts were 
extracted using two circular regions of radius 15" located 
in source-free regions. The count rates were obtained 
from the image lists using the Swift tool uvotsource. 
As GRB 2210094 fades, even using a small aperture, 
the contribution from the nearby star increases in dom- 
inance. In order to estimate the level of contamination, 
for each filter we combined the late time exposures be- 
tween To + (3.4 x 109) s and To + (4.4 x 109) s. We ex- 
tracted the count rate in the late combined exposures 
using the same 2.5" aperture and applied an aperture 
correction to 5". These were subtracted from the source 
count rates to obtain the afterglow count rates. The af- 
terglow count rates were converted to magnitudes using 
the UVOT photometric zero points (Poole et al. 2008; 
Breeveld et al. 2011), as well as to fluxes. To improve 
the signal-to-noise ratio, the count rates in each filter 


were binned using At/t —0.2. The final photometry is 
plotted in the bottom panel of Figure 2 and listed in 
Appendix C. 


2.4. MAXI 


MAXI, which is onboard the International Space Sta- 
tion (ISS), performs scanning observations of about 8096 
of the sky every 92 min. Prior to the GBM detection of 
the prompt emission from GRB 221009A, the two Gas 
Slit Camera (GSC; Mihara et al. 2011) units, GSC. 4 and 
GSC. 5, scanned the source region at 12:25 UT on 2022 
October 9 (15 —3.1 ks), and no enhancement was seen at 
that time. The 2-10 keV 30 upper limit is 0.074 ct cm~? 
s71, corresponding to a flux of approximately 7.9x 10-1? 
erg em ? s Í. At 13:58 UT (To4-2.5ks) the GSCs de- 
tected a very bright transient with a 2-10keV flux of 
5.97 + 0.19 ct cm? s^! (~ 6.5x10 8 erg cm? s^ 5; 
Negoro et al. 2022). Unfortunately, at the time the 
data were collected, there was a loss of the signal to 
the ground, and the data were stored in the High Rate 
Communication Outage Recorder in the ISS, leading to 
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a delay in reporting the transient. Before downloading 
these stored data at around 15:55 UT, the MAXI/GSC 
Nova-Alert System (Negoro et al. 2016) triggered on the 
source at 15:31:03 UT (7o--8.0 ks), about 81 min after 
the BAT trigger. In this second observation by MAXI, 
the 2-10 keV flux had dropped to 1.00£0.08 ct em ? s71 
(about 1x10 5 erg cm"? s^ 1). 

Thanks to the exceptional brightness, the GSCs sig- 
nificantly detected the afterglow during four scans that 
followed the second detection, to 21:41 UT, six scans in 
total (e.g., Kobayashi et al. 2022). Typically, GRBs are 
not detected beyond the first one or two scans, taken at 
~ 92 min intervals . After the six detections, the sky re- 
gion containing the source was not observable by MAXI 
until 15:31 UT on 2022 October 11. 

We performed a spectral analysis using data obtained 
from GSC.4 and GSC.5 (for details see Appendix D). 
We first attempted a standard analysis using the same 
model applied to the XRT spectra, i.e., a power-law with 
two absorption components, with a fixed Galactic ab- 
sorption column density of Ng = 5.38 x 10?! em ?. 
However, we obtained steeper photon indices and/or 
lower intrinsic absorption column densities than those 
of the XRT (see Table 6). If we fixed Nm;aps at 
the weighted mean value obtained with the XRT at 
To-- (3.4 — 50) ks, 1.29x 107? cm ?, even steeper photon 
indices were obtained. 

'These discrepancies can be explained by the GSC 
spectra containing dust scattered soft X-rays, which can 
not be spatially resolved by the GSC. Mitigation of the 
dust scattering process is discussed in Appendix D; to 
summarize, an additional spectral component to account 
for the scattered emission was included with a power-law 
index derived from XRT PC-mode data and the associ- 
ated flux allowed to vary as a free parameter. The af- 
terglow photon index was fit as a free parameter in the 
data from the first two scans, but fixed at the weighted 
mean XRT value at To+(10 — 50) ks due to poor statis- 
tics. The resulting spectral parameters are displayed 
in Table 1, and flux measurements (extrapolated to the 
0.3-10.0 keV band) are plotted in Figure 2. 


2.5. NICER 


NICER received notification of this GRB from MAXI 
on the ISS via OHMAN (On-orbit Hookup of MAXI and 
NICER) at 14:10:57 UT, but could not begin regular 
observations due to poor visibility until Tọ + 52.870 ks. 
Only two observations were possible during this early 
phase. The first and second observations were at 7o + 
14.003 ks and Tọ + 19.560 ks with an exposure time of 
29 s and 62 s, respectively. After To + 52.870 ks, NICER 


intermittently made observations over a period of about 
42 days. 

NICER’s X-ray Timing Instrument (XTI) consists of 
56 modules, each comprising an X-ray concentrator op- 
tic and silcon drift detector (Gendreau et al. 2012). Dur- 
ing these observations 52 X' TI modules were operational. 
We reduced the data with the nicerl2 command and 
generated level 2 cleaned events. To avoid potential con- 
tamination of instrumental noise, we further excluded 
XTI module numbers 14 and 34. We defined a time in- 
terval in which an exposure time of more than 200 s can 
be continuously secured as a one-time interval and cre- 
ated an energy spectrum at each time interval (except 
for the first two observations) with good statistics. 

'To select a time period with less background influence 
due to charged particles, we also checked the distribu- 
tion of the rate of 12-15 keV photons, where the op- 
tics have less effective area. Most data are distributed 
around 0.048 ct s^! and can be represented by a nor- 
mal distribution with a standard deviation of 0.018 ct 
s Í. Therefore, we removed data with a rate greater 
than 20 from the mean in the 12-15 keV range after 
794-100 ks seconds where the background contribution 
becomes large (> 0.084 ct s Í). We also did not use 
data after 794-1000 ks because the background rate in- 
creased due to a geomagnetic storm. The background 
spectral model is estimated by the 3C50 background 
model (Remillard et al. 2022) using the nibackgen3c50 
command. 

'To explore the spectral evolution as in 82.2, we di- 
vided the NICER data into six time intervals and fitted 
the time-resolved spectra with a power-law model and 
the same two absorption components. Since NICER is 
a non-imaging detector, it cannot remove photons orig- 
inating from the dust echo within its 5’ diameter field- 
of-view. To evaluate the impact of the dust echo, we 
fitted the data in two energy ranges, 1-8 keV and 4- 
8 keV bands. Since the data fitting in the 4-8 keV range 
does not allow us to constrain the intrinsic Ng value, we 
fixed it to 1.29x 10?? cm^? which is the weighted mean 
value obtained by XRT at 79--(10.0—961.6) ks (Table 1). 
Comparing the fitting results in the two energy bands, 
the results for the 4-8 keV band show that the photon 
indices are systematically harder than those for the 1- 
8 keV band results and consistent with the XRT results. 
Therefore, we conclude that 4-8 keV light curve is much 
less affected by the dust echo, and report these results 
in Table 1. 


3. ANALYSIS 
3.1. Dust Scattering Echo 


While complicating the point source analysis of the 
afterglow emission, the dust scattering echo contains a 
rich amount of information about the intervening Galac- 
tic dust along the line of sight to GRB 2210094. Dust 
echoes from GRBs offer the unique opportunity to mea- 
sure distances to dust clouds geometrically, with a pre- 
cision limited only by the angular resolution of the tele- 
scope. This is because the distance Darg to the X-ray 
source is effectively infinite, relative to the distance to 
the dust. This allows us to eliminate the source distance 
from the scattering angle equation (e.g., Heinz et al. 
2015) and solve it directly for the distance Daust to the 
dust (using the small angle approximation): 


Dana 2cAt 
T (1) 


where At is the time delay between the GRB and the 
exposure and @ is the off-axis angle of the ring relative 
to the position of the GRB. 

Because the vast majority of the fluence of the GRB 
was concentrated in the prompt GRB emission, to lowest 
order, we can approximate the input light curve that 
caused the echo as a delta function in time. This implies 
that there is a single, unique At for each photon detected 
by the XRT, and thus, each photon can be referenced to 
a specific dust distance in eq. 1. This makes the analysis 
of this echo substantially simpler than in comparable 
echoes from Galactic X-ray transients. 

Details of our analysis of the dust echoes from 
GRB 221009A are provided in Appendix E. Figure 4 
shows the best fit column density histogram for a 
standard Mathis, Rump, & Nordsieck (MRN) model 
(Mathis et al. 1977), assuming a soft X-ray fluence of 
Fo.s-5kev = 2.1 x 107? erg cm7? (Lesage et al. in prep.). 
The figure shows clear evidence for several dust compo- 
nents on the near side of the Galaxy, with the largest 
column densities found within distances of about 1 kpc, 
as expected given the Galactic latitude of GRB 2210094. 

However, there is also clear evidence for dust at larger 
distances, located well above the Galactic disk (see 
also Negro et al. 2023). While the column densities of 
these dust concentrations are small, the large cross sec- 
tion at the correspondingly small scattering angles of 
this more distant dust makes the scattering from these 
clouds detectable in the radial intensity profiles, as can 
be seen in Figure 13, for example. As such, the dust 
echo from GRB 221009A proves that echo tomography 
is exquisitely sensitive to small dust concentrations at 
large distances that are not typically measurable with 
other techniques. 

'The largest column density cloud is located approx- 
imately between 400 and 600 pc distance, with a very 
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Figure 4. Dust column density (expressed as equivalent 
hydrogen column density, assuming solar abundances) de- 
termined from fitting dust scattering models to the radial 
intensity profiles as a function of distance, for a standard 
MRN model. The dependence on the uncertain soft X-ray 
fluence of GRB 221009A is shown in the y-axis label. Vertical 
gray lines indicate lo error bars. 


nearby dust cloud at a distance of about 200 pc and 
a third cloud at distance approximately 700 pc. Nat- 
urally, these structure are very local compared to the 
large scale Galactic structure. The very rapid detection 
of the echo by the XRT allows us to map very nearby 
dust typically hard to study using dust echoes. 

We find a total Galactic column density of Ng ~ 5.9 x 
10?! em ?, broadly consistent with the Galactic column 
density of 5.38 x 10?! em ? found by Willingale et al. 
(2013). Given that the largest uncertainty in this value 
derives from the poorly constrained estimate of the soft 
X-ray fluence of GRB 2210094, we can infer that our 
assumed fluence of Fo.s—5key = 2.1 x 107? ergscm ? is 
likely correct within roughly 10%. 

We find that a simple MRN distribution provides an 
excellent fit to the temporal evolution of the radial in- 
tensity profile, without the need for more complex dust 
chemistry (such as ice mantles). In this sense, the dust 
towards GRB 221009A appears to be typical for inter- 
stellar dust seen along sightlines to other recent X-ray 
light echoes, such as the echo from V404 Cyg (Beard- 
more et al. 2016; Heinz et al. 2016; Vasilopoulos & 
Petropoulou 2016). 


3.2. Broadband Spectral Energy Distribution 


The combined XRT+MAXI+NICER dataset allows 
us to explore the soft X-ray afterglow evolution in great 
detail. In the earliest soft X-ray spectra (To 4-2.5-10 ks), 
the afterglow exhibits strong spectral evolution, with 
the XRT WT-mode power-law photon index increasing 
(i.e., getting softer) by = 1796 (Table 1). This change 
is much larger than the inferred WT-mode systematic 
uncertainty (Appendix B), and is corroborated by the 
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2-20keV MAXI spectra derived around this time (Ap- 
pendix D). 

From To + 10-50 ks, the spectra show no definitive 
signs of evolution; there are hints at a reduction in the 
host galaxy absorption throughout the XRT W'T-mode 
data, although this trend is reversed when XRT switched 
to PC mode. Given the differences between the WT- and 
PC-mode results are much less than the dust-induced 
systematics on the WT results (Appendix B), and no 
similar spectral evolution is observed in the NICER 
spectra, we do not believe there is supportable evidence 
for spectral evolution at these times. 

On the longer timescales probed by the XRT PC-mode 
data, clear spectral evolution is seen and takes the form 
of an increasing (i.e., softer) photon index, which is con- 
firmed by the 4-8 keV NICER spectra (Table 1). 

To constrain the broadband spectral behavior, we 
built and fit two spectral energy distributions (SEDs) 
at To + 4.2ks and To + 43ks, following the procedure 
outlined in Schady et al. (2010, 2007), using data from 
the BAT survey mode, XRT and UVOT. For the 4.2 ks 
SED we use data in the range 2.5-4.7 ks. For the 43ks 
SED we used data in the range 23-103 ks. Details of the 
SED construction are provided in Appendix F. 

'The SEDs were fit using XSPEC. We tested two differ- 
ent models for the continuum: a single power-law, and a 
broken power-law with the change in spectral slope fixed 
to be AG = 0.5 (corresponding to the expected change 
in spectral slope caused by the synchrotron cooling fre- 
quency; Sari et al. 1998). In each of these models, we 
also included two dust and gas multiplicative compo- 
nents to account for the Milky Way and host galaxy 
dust extinction and photoelectric absorption (TBABS, 
ZTBABS, ZDUST)*. We also include a component for 
attenuation by the intergalactic medium (zIGM). The 
Galactic components were frozen to the reddening and 
column density values from Schlegel et al. (1998) and 
Willingale et al. (2013), respectively. In the analysis, 
we found consistent results between Milky Way, Small 
Magellanic Cloud (SMC) and Large Magellanic Cloud 
dust extinction laws, therefore we only report the fits 
using the SMC dust extinction law. 

Due to the uncertain contribution to the XRT WT 
spectra by the dust rings, we extracted W'T spectra us- 
ing the 5-20 pixel annular extraction region identified in 
Appendix B for the spectrum for the early SED, and a 10 
pixel radius circular extraction region for the later time 
SED, chosen to minimize the fraction of dust scattered 


^ We used ZDUST for both the Galaxy and host galaxy dust com- 
ponents, but with the redshift set to zero for the Milky Way dust 
component. 


emission in the extraction region, with a background 
spectrum taken from a 75-125 pixel annular region.? 
Since the intrinsic Ny in the WT fits is sensitive to the 
dust contamination (see Appendix B and Table 5), we 
also provide fits using BAT and UVOT alone. The re- 
sults of the analysis are given in Table 2 and plotted in 
Figure 5. 

Examining the SED at 754-4.2 ks including BAT, XRT 
and UVOT data, we find that a broken power-law is pre- 
ferred, with the F-test finding a statistical improvement 
of > 50. The break energy required is at = 7 keV. At 
the same epoch, excluding the XRT data we find again 
a broken power-law model is preferred, with the F-test 
suggesting the improvement is > 3c. The photon index 
for the broken power-law model is consistent with the fits 
including the XRT spectral files. However, the intrinsic 
dust extinction is lower: 0.23 + 0.05 mag compared to 
0.51 + 0.03 mag including all data. Not surprisingly, ex- 
cluding the XRT data greatly increases the uncertainty 
in the break frequency. 

For the SED at To +43 ks, a broken power-law is again 
preferred, with the F-test suggesting a break is statisti- 
cally required at > 30. The break frequency remains in 
the X-ray bandpass, with evidence for a decline in energy 
over this period. The fits still prefer a significant host 
galaxy reddening, but the exact value of E(B — V)host 
again varies depending on the inclusion or exclusion of 
the XRT data. 

To summarize, both at To + 4.2 ks and To + 43 ks, we 
find strong evidence favoring a broken power-law model, 
with the lower energy photon index T ~ 1.7. Both 
epochs favor a break energy in the X-ray bandpass, and a 
modest amount of intrinsic absorption/reddening in the 
host galaxy: E(B — V)host © 0.3-0.5 mag; Ng ~ (1.1- 
1.4)x10?? em ?. More precise constraints, however, are 
precluded due to the complications resulting from the 
dust scattering echoes. 


3.3. Light Curve 


The combined X-ray and UV/optical light curves of 
GRB 221009A are shown in Figure 2 (top and bottom 
panels, respectively.. The XRT, MAXI, and NICER 
measurements are placed in a common bandpass from 
0.3-10.0 keV. The BAT survey mode data, however, are 
not extrapolated to this band and instead plotted from 
14-195 keV — due to the inference of a spectral break 


5 See https:/ /www.swift.ac.uk/analysis/xrt/backscal.php. 
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Table 2. Fits to the two spectral energy distributions, at To + 4.2ks and To + 43ks, built using UV/optical (U), X-ray (X) 
and gamma-ray (B) data. Both epochs were fit with a power-law (POW) and a broken power-law (BKPOW) continuum, 
accounting for Galactic and host galaxy gas and dust. The columns are time of the SED; model; host galaxy equivalent 
column density Nz,x,host; host galaxy reddening, E(B-V)nosc; power-law photon index or the first photon index of the broken 
power-law model I, (T2 is fixed to be [ + 0.5); break energy for the broken power-law spectral models Ex, x? & degrees of 


freedom (d.o.f); and the null hypothesis probability. 


SED Model N n, x host E(B-V)nost T Eor x? (d.o.f) | Null Hypothesis 
(10 cm~?) (mag) (keV) Probability 
B+U 4.2 ks POW = 0.26+0.05 1.74+0.02 = 16 (9) 7.3e-02 
B+U 4.2ks | BKPOW - 0.234-0.05 — 1.707003 dig. 4 (8) 8.4e-01 
B+X+U 4.2ks POW 1.75 +0.03 1.04 0.02 1.93 0.01 - 1372 (681) 6.9e-49 
B+X+U 4.2ks BKPOW 1.35+Æ0.03 0.51+0.03 1.69 0.01 6.8 0.3 861 (680) 2.7e-06 
B+U 43ks POW = 0.38 +0.06 1.89 + 0.04 = 14 (9) 1.2e-1 
SS) 

B+U 43 ks BKPOW = 0.38 +0.06 1.39+0.04 73x IT ipis ue 14 (8) 8.3e-02 
B--X--U 43 ks POW 1.27 +0.03 0.492: 0.03 1.83 0.01 = 778 (678) 4.2e-03 
B+X+U 43ks BKPOW 1.1440.08 0.28+0.03 1.72 4 0.02 4.5+ 0.2 748 (677) 3.0e-02 
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Figure 5. Spectral energy distributions at To  4.2 ks and To 43 ks, containing BAT, XRT and UVOT data. The top panel 
displays the best fit unabsorbed broken power-law model (red) together with the data (black) and the data folded with the 
model (orange). The bottom panel shows the ratio between the data and the model folded through the instrument response. 


below this band (83.2), these points are instead meant 
to be illustrative. 

While soft X-ray observations of the afterglow did not 
begin until To +2.5 ks, the light curve clearly lacks signa- 
tures of the canonical X-ray afterglow behavior at early 
times, including the prompt decay phase, the plateau 
phase, and any prominent flaring (Nousek et al. 2006). 
'This relatively simple monotonic decay has been noted 
previously in other GRBs detected at high (>GeV) en- 
ergies (Yamazaki et al. 2020). 


We fit the joint 0.3-10.0 keV (XRT+MAXI+NICER) 
X-ray light curve using the Markov Chain Monte Carlo 
software package EMCEE (Foreman-Mackey et al. 2019). 
An additional fixed fractional error was included as 
a nuisance parameter to account for cross-calibration 
uncertainty. We find a broken power-law model is 
strongly preferred over a single power-law, although still 
formally a more complex behavior is indicated. The 
best-fit parameters for the broken power-law model are 
Qi, = —1.498 + 0.004, ag x = —1.672 + 0.008, and 
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tbreak,x = (7.9* 15) x 104s (68% confidence intervals). 
The initial MAXI detection at To  2.5ks clearly indi- 
cates a shallower decay at early times (To4- < 3.3 ks), 
while a late time excess (To + 2 Ms) may indicate some 
energy injection (e.g., Zhang et al. 2006). 

For the UV/optical data, we constructed a combined 
single-band light curve following the procedure outlined 
in Schady et al. (2010). We then fit this light curve using 
the identical procedure described above for the X-ray 
data (including the cross-calibration nuisance parame- 
ter). We find that a single power-law and broken power- 
law model provide comparable quality fits to the data — 
the statistically preferred model depends sensitively on 
how late-time (> 10? ks) measurements are included in 
the fit. The best fit index derived for the single power- 
law model is ag = —1.13 + 0.01, while for the broken 
power-law model the parameters are ao; = —0.98*0 bs, 
ios = LAICUS. and (aso = O27) x 104s. 
Regardless of the model selected, the UV /optical light 
curve clearly declines more slowly than the X-ray emis- 
sion. Given the large foreground extinction, these wave- 
lengths are unaffected by contamination from an emerg- 
ing supernova or an underlying host galaxy, and thus 
this accurately reflects the afterglow behavior. 


4. DISCUSSION 
4.1. Comparison with previous GRBs 


The top panel of Figure 6 plots the observed 0.3- 
10.0 keV flux light curve of GRB 221009A, overplotted 
on the entire sample of XRT light curves since Swift 
launch (all uncorrected for absorption). At the time of 
the first XRT observations, GRB 221009A is approxi- 
mately an order of magnitude brighter than any previ- 
ous X-ray afterglow. Put differently, forming the dis- 
tribution of X-ray afterglow flux at To + 4.5 ks (a time 
selected to minimize the number of events requiring in- 
terpolation due to orbital gaps in light curve coverage), 
GRB 2210094 falls 5.50 above the mean. 

However, GRB 221009A is also one of the nearest 
GRBs detected since the launch of Swift — 12 events 
(out of 22400) have lower redshifts reported in the online 
Swift GRB Table?. To compare to the intrinsic prop- 
erties of the X-ray afterglow sample, we converted the 
observed count-rate light curve into intrinsic luminosity 
in the 0.3-10 keV energy band in the GRB comoving 
frame. The luminosity in a given bin is given by: 


L = 4r DŻ RC,k (2) 


6 See https://swift.gsfc.nasa.gov/archive/grb table. 


where Dz is the luminosity distance (749.3 Mpc), R 
is the measured 0.3-10 keV count rate (in the observer 
frame), C, is the conversion from count-rate to unab- 
sorbed 0.3-10 keV flux (also in the observer frame) ob- 
tained from the spectral fits above (82.2), and k is the 
k-correction of Bloom et al. (2001), which corrects from 
observed 0.3-10 keV flux to the 0.3-10 keV flux in the 
GRB comoving frame. To calculate k, we assume the 
spectrum is an unbroken power-law with the photon in- 
dex l = 1.8. The time axis must also be corrected to the 
observer frame, by dividing the time by 1+ z = 1.151. 
'The resultant light curve is shown in the bottom panel 
of Figure 6. 

While the lack of early time observations some- 
what complicates the comparison, it is clear that 
GRB 221009A is at the high end of the X-ray after- 
glow luminosity distribution (1.80 above the mean). But 
unlike when considering this comparison in flux space, 
GRB 221009A is by no means an outlier. Rather we 
conclude that the exceptional observed brightness of the 
X-ray afterglow of GRB 2210094 results from the com- 
bination of being very nearby and very intrinsically lu- 
minous: while neither the distance nor the luminosity 
are unprecedented, together they make GRB 221009A 
unique amongst the Swift afterglow sample. 

Comparing GRB221009A to the broader popula- 
tion of UV/optical afterglows is further complicated 
by the large (and uncertain) line-of-sight extinction, 
both foreground and in the host galaxy. With my = 
17.59*013 mag at To + 3574s, GRB 2210094 is signifi- 
cantly fainter than the brightest UVOT-detected after- 
glows, such as GRB 080319B (my ~ 15.4 mag; Racusin 
et al. 2009; Bloom et al. 2009) and GRB 130427A (mu = 
13.9 mag; Maselli et al. 2014) at comparable observer- 
frame times post-burst. But even applying only a correc- 
tion for foreground (i.e., Milky Way) extinction of Ay = 
6.8mag (Schlafly & Finkbeiner 2011), the UV/optical 
emission from GRB 221009A would likely have been the 
brightest observed to date had it occurred at high Galac- 
tic latitude. 

If we adopt E(B — V)host = 0.4mag and an SMC- 
like extinction law (83.2), we infer a U-band host ex- 
tinction of 1.8 mag. Neglecting k-corrections (which are 
likely to be much smaller than the uncertainty in the 
extinction correction), we find an absolute magnitude 
of My ~ —29.3mag (AB). Comparing to Figure 7 in 
Perley et al. (2014), we find that the UV /optical light 
curve of GRB 2210094 occupies a similar location in lu- 
minosity space as the X-ray afterglow: towards the most 
luminous end, but entirely consistent with the existing 
distribution. 
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Figure 6. A comparison of GRB 2210094 with all of those 
observed by XRT. The greyscales indicate the number of 
GRBs in each (time, brightness) bin, the blue light curve is 
GRB 2210094, with two notable bright GRBs, GRB 130427A 
(purple) and GRB 080319B (orange) shown for comparison. 
Top: Observer frame, 0.3-10 keV observed flux light curves. 
Bottom: Comoving-frame, 0.3-10 keV intrinsic (unabsorbed) 
luminosity light curves; the greyscale sample data includes 
only those GRBs with published redshifts. 


4.2. Astrophysical Rate 


Although Swift was behind the Earth at To, we use the 
GBM light curve and spectrum (Lesage et al. in prep.) 
to estimate what GRB 221009A would have looked like 
to the BAT, as well as to determine the occurrence rate 
of such events. Adopting the time-averaged spectrum 
of the main pulse from To + 218.5-277.9s (10-1000 keV 
fluence of 8.293 x 107? erg em ?), we derive a BAT 
(15-150 keV band) fluence of 1.86 x 107? erg em ? s^! 
and an average flux of 3.13 x 107^ erg em ? s-!. This 
fluence from the main pulse alone is 50x higher than the 
largest fluence detected to date by BAT (GRB 1304274; 
Maselli et al. 2014). 
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At z — 0.151, the 15-150 keV (rest-frame) isotropic 
prompt energy release for the main pulse is Fy iso = 
9 x 105? erg. We employ the 15-150 keV energy range to 
facilitate a direct comparison with the entire population 
of BAT GRBs with measured redshifts. 

'To estimate the relative occurrence rate of such an en- 
ergetic event, we compare GRB 2210094 with the intrin- 
sic long GRB luminosity function derived in Lien et al. 
(2014). We consider only the main pulse, as this domi- 
nates the burst energetics, and it is difficult to ascertain 
when the prompt emission ends, given the afterglow was 
bright enough to trigger the BAT at To + 3.4 ks. 

We randomly generated 104 GRBs using the intrin- 
sic GRB rate and luminosity distribution in Lien et al. 
(2014), and each simulated burst was randomly assigned 
a pulse structure drawn from the real BAT GRB sam- 
ple. We calculated the E} iso of these simulated bursts, 
and found that only one burst had an energy release 
(slightly) higher than the main pulse of GRB 221009A. 
Therefore, we conclude only ~1/10* long GRBs are as 
energetic as GRB 221009A. 

Using this value, we crudely estimate the rate at which 
such energetic GRBs are detectable by BAT: 


NBAT RI Hana fEiso frov Ísurvey = 0.06 yt (3) 


where Rears is the intrinsic all-sky long-GRB rate from 
Lien et al. (2014), fgio = 1.0 x 1074 is the fraction 
of GRBs in this intrinsic sample that have Fy iso larger 
than GRB 221009A, froy = 1/6 is the fraction of sky 
covered by the BAT field of view, and fsurvey = 0.8 is 
the fraction of time BAT is able to trigger on GRBs (ne- 
glecting spacecraft slews, SAA passage, etc.). In other 
words, we need to wait = 1/0.06 = 17 yr for an event 
like the main pulse of GRB 2210094 to occur in the BAT 
field of view. 

We emphasize, however, that the above estimate sim- 
ply indicates that BAT should have detected ~1 GRB 
more energetic than GRB 2210094 over its lifetime, in- 
dependent of distance. But in addition to being highly 
energetic, GRB 2210094 is also one of the most nearby 
GRBs detected by Swift. To estimate the intrinsic (vol- 
umetric) rate of comparable events, we utilize the BAT 
trigger simulator (Lien et al. 2014) to calculate the de- 
tectability of the prompt emission at different redshifts 
under different representative geometries. The full de- 
tails of the simulation are provided in Appendix G. 

Using this framework, we derive an upper limit on 
the local volumetric rate of GRB 221009A-like events of 
RanB,comov(2 = 0) < 6.1 x 107* Gpc^? yr 1. Integrat- 
ing this flat comoving rate from z — 0 to z — 10, we 
obtain an upper-limit on the all-sky rate of such events 
of < 0.5yr-!. Comparing to the all-sky intrinsic long 
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GRB rate of ~ 4571 yr Í in Lien et al. (2014), the frac- 
tion of GRB 221009A-like events is roughly 0.5/4571 < 
1.0 x 1074, which is similar to the relative rate derived 
above. 

Even more remarkable, we can use these results to 
derive the rate of GRB 221009A-like GRBs within the 
volume out to z = 0.151 (1.1 Gpc ?). This implies 
we would need to wait over 210? years to detect 
another GRB 221009 A-like event within this vol- 
ume. The combination of the large energy release and 
the small distance make GRB 2210094 truly a once in a 
lifetime phenomenon. 


4.3. The Nature of Energetic GRBs 


Interpreting GRB 2210094 in the context of the stan- 
dard afterglow synchrotron model (e.g., Sari et al. 1998) 
presents a number of challenges, in particular if the 
X-ray and optical emission is assumed to arise from a 
common origin (c.f., Ghisellini et al. 2007; de Pasquale 
et al. 2009). To begin with, the broadband spectral 
fitting performed in 83.2 strongly favors the presence 
of a spectral break around the soft X-ray band, with 
the change in spectral slope being consistent with the 
cooling frequency ve. However, if we adopt relatively 
standard afterglow parameters (eg = 0.01, p = 2.5, 
Ny = 0.15), then we find that a cooling frequency = 
5 keV (1015 Hz) at ~0.1d post-trigger implies a jet ex- 
panding into an extremely low density circumburst envi- 
ronment: ng = 107? cm ^? for a constant-density (ISM- 
like) environment. Such low densities have been inferred 
previously for highly energetic events (e.g., Cenko et al. 
2011), but this remains difficult to reconcile with the 
massive star progenitors of long GRBs. 

Furthermore, standard afterglow closure relations 
(e.g., Racusin et al. 2009) are unable to reproduce the 
observed spectral and temporal power-law indices. For 
example, the broadband SED fits presented in 83.2 find 
an optical spectral index fo æ 0.7. This is consistent 
with the temporal decay observed for a jet expansion 
into a constant-density medium (ag ~% 1.1; pæ 2.4). 
However, this would predict a significantly shallower 
temporal decay in the X-rays (ax ~ 1.3) than what 
is observed. While expansion into a wind-like medium 
would predict a steeper X-ray temporal decay, the opti- 
cal emission is inconsistent with this picture. 

Given the exceptionally large Ey iso, GRB221009A 
requires a narrow opening angle (and hence an early 
jet break time) in order to be compatible with the 
geometry-corrected energy release inferred from the pre- 
Swift sample (e.g, Frail et al. 2001). Assuming an effi- 
ciency of x15% in converting bulk kinetic energy into 
prompt y-ray emission, a beaming-corrected energy re- 


lease of < 10°” erg requires an opening angle of 0; < 4°. 
For a top-hat jet expanding into a constant-density 
medium, this implies a jet break time of: 


ua (4) 


t; < 3. — 
j $355 E 


where no is the circumburst density in units of 0.1 cm7? 
(Sari et al. 1999). Note that due to the -1/3 exponent, 
this result is relatively insensitive to the exact value of 
the circumburst density. 

With X-ray coverage extending out to To + 70d, the 
Swift+MAXI+NICER X-ray light curve is more than 
adequate to search for such a geometric signature. But 
out to late times, the temporal decay does not steepen 
beyond ax = 1.7, significantly shallower than the a > 2 
behavior anticipated for a simple on-axis top-hat jet. 

The X-ray afterglow steepening observed at To + 8 x 
10* s could conceivably be attributed to a jet break (see 
also D'Avanzo et al. 2022). In this case, the corre- 
sponding jet would indeed be extremely narrow, with an 
implied opening angle of ~ 2°. However, as described 
above, the post-break decay index is significantly shal- 
lower than expected for an on-axis top-hat jet. Interest- 
ingly, similar behavior was observed in GRB 130427A 
and attributed to time-variable microphysical parame- 
ters (Maselli et al. 2014). One alternative explanation is 
a complex angular structure for the jet (i.e., compared 
to the simple top-hat model assumed above, with a con- 
stant energy as a function of angle). The presence of 
energetic material outside the jet “edges” would nat- 
urally account for the shallower decline initially (e.g., 
O'Connor et al, in prep.). However, eventually when 
viewed sufficiently far off-axis the decay should steepen, 
something not observed in the X-ray light curve out 
to To + 70d. A detailed exploration of the jet struc- 
ture would require additional multi-wavelength obser- 
vations, in particular radio coverage, and is thus beyond 
the scope of this work. 


5. CONCLUSIONS 


Here we present the combined Swift, MAXI, and 
NICER observations of GRB 2210094, spanning from To 
+ 2.5 ks to To + 73 days and from UV to hard X-ray en- 
ergies, providing the opportunity to study the prolonged 
afterglow in exquisite detail. The primary conclusions 
are as follows: 


e Through dust echo tomography, we identify a se- 
ries of dust clouds along the line-of-sight in our 
galaxy at distances ranging from 200 pc to beyond 
10 kpc (i.e., well above the Galactic disk). As such, 
the dust echo from GRB 2210094 proves that echo 


tomography is acutely sensitive to small dust con- 
centrations at large distances that are not easily 
measurable with other techniques. 


e The X-ray afterglow of GRB 221009A is more than 
an order of magnitude brighter at To + 4.5 ks (in 
terms of measured flux) than any previous GRB 
observed by Swift. However, the intrinsic X-ray 
afterglow luminosity is at the high end, but con- 
sistent with, the distribution of Swift GRBs with 
redshifts. Thus, GRB 2210094 is unique both be- 
cause it is luminous and very nearby (by GRB 
standards). 


e We calculate the rate of GRB 221009A-like events 
(ie, GRBs with prompt emission as ener- 
getic/luminous as that inferred by the GBM). We 
find that only =1 in 10^ GRBs has an Ey iso com- 
parable to GRB 221009A. When factoring in the 
distance, we find that the rate of GRBs as ener- 
getic and nearby as GRB 221009A is < 1 per 1000 
years. 


e From the extensive multi-wavelength data sets, 
we find that the afterglow emission from 
GRB 221009A is not well described by standard 
synchrotron afterglow theory. In particular, the 
break observed in the X-ray light curve at To + 
79 ks is inconsistent with a jet break from an on- 
axis top-hat jet. Either the jet structure is more 
complex, with significant power outside the jet 
core, or GRB 2210094 is not narrowly collimated. 
'The latter scenario would have profound implica- 
tions for the energy budget of the event. 


While GRB 221009A disappeared behind the Sun 
for most observatories at ~70d post-burst, the X-ray 
and radio afterglow are still anticipated to be eas- 
ily detectable after the field becomes visible again in 
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February 2023. Coupled with the exquisite multi- 
wavelength data collected from radio to VHE gamma- 
rays, the broadband story of this exceptional event will 
continue to unfold over the coming months, years, and 
(perhaps) decades. 
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APPENDIX 


A. BAT SURVEY MODE AND MOSAIC FLUX MEASUREMENTS 


Details of the flux measurements and upper limits obtained from the Swift BAT survey mode data are shown in 
Table 3, and the results of the daily mosaicing are presented in Table 4. 


B. Swift XRT DUST AND BACKGROUND SUBTRACTION 


In WT mode the Swift XRT CCD is read out in columns and the spatial information collapsed to one single dimension. 
As a result, the ‘source’ region will contain all dust photons ‘above’ and ‘below’ the GRB in terms of CCD position, 
while the background regions will sample a different vertical slice through the dust rings, with potentially different 
spectral properties from those in the source region. Furthermore, the columns are summed over the full 600-pixel 


Table 3. BAT Survey Observations of GRB 221009A. Upper limits are calculated with an assumed photon index of 1. 


Obs ID Pointing ID Tstart — To UTC Time Exposure 14-195keV Flux Photon Index 
(ks) (ks) (erg cm? s^!) 
03111868007 20222810634  —110.5 2022 October 08 06:34:32 0.661 « 1.11 x 1075 
20222821110 -7.6 2022 October 09 11:10:19 0.373 < 1.55 x 1075 — 
00015314118 20222820117 — —43.2 2022 October 09 01:17:08 0.924 < 6.33 x 107? — 
00046390010 20222820135 — —42.1 2022 October 09 01:34:40 0.172 < 3.34 x 1075 — 
00015314119 20222820626 ^ —24.6 2022 October 09 06:26:19 0.900 < 6.36 x 107? — 
00015357004 20222820921 ^ —14.1 2022 October 09 09:21:28 ^ 0.544 < 1.42 x 1075 — 
00011105067 20222821046 —94 2022 October 09 10:46:08 1.344 < 1.36 x 1075 — 
00015314120 20222821248 -1.7 2022 October 09 12:48:26 0.846 < 7.14 x 107° — 
01126854000 20222821422 3.9 2022 October 09 14:22:00 0.627 . 3.227032 x 1078 2.131943 
01126853001 20222821917 21.6 2022 October 09 19:17:08 0.500 — 2.61*092 x 10-? 1.811944 
01126853003 20222822027 25.8 2022 October 09 20:26:45 0.900 — 2.44*069 x 107° FL í 
00015314121 20222822047 27.0 2022 October 09 20:46:57 . 0.755 . 3.01*14i x 107° 1.587021 
01126853004 20222822156 31.1 2022 October 09 21:55:37 1.744 . 1.73*0$) x 107° 1.854939 
20222822332 36.9 2022 October 09 23:31:40 1.500 < 4.09 x 107° — 
20222830129 43.9 2022 October 10 01:29:08 1.073 < 4.76 x 107? — 
20222830243 48.4 2022 October 10 02:43:05 1.742 < 3.82 x 107° — 
20222830417 54.0 2022 October 10 04:16:45 1.738 < 3.74 x 107? — 
20222830754 67.0 2022 October 10 07:53:52 0.914 < 5.05 x 107? — 
20222830934 73.0 2022 October 10 09:34:21 0.514 < 6.78 x 107° — 
20222831234 83.8 2022 October 10 12:33:38 1.314 < 4.61 x 107? — 
20222831354 88.6 2022 October 10 13:53:49 1.183 < 4.81 x 107? = 
20222831531 94.4 2022 October 10 15:30:56 ^ 0.900 < 5.34 x 107? — 
20222831550 95.6 2022 October 10 15:50:56 0.543 < 6.53 x 107? — 
20222831704 100.0 2022 October 10 17:04:26 1.679 < 3.76 x 107° — 
03111808004 20222830447 55.8 2022 October 10 04:46:38 ^ 0.594 < 5.99 x 107° — 
00033349157 20222830906 71.3 2022 October 10 09:05:56 ^ 1.655 < 3.68 x 107° — 
00015314124 20222831417 90.0 2022 October 10 14:17:05 0.807 < 7.28 x 107° — 
01126853005 20222831837 105.6 2022 October 10 18:37:15 1.500 < 3.80 x 107? — 
00015314125 20222831904 107.2 2022 October 10 19:03:48 0.704 2711 x 107° = 
01126853006 20222832022 111.9 2022 October 10 20:22:19 1.500 « 4.07 x 107? — 
20222832147 117.0 2022 October 10 21:47:19 0.900 < 5.17 x 107? — 
20222832204 118.2 2022 October 10 22:07:19 0.299 < 8.71 x 107? — 
20222840108 129.1 2022 October 11 01:08:33 0.839 < 5.31 x 107? — 
20222840233 134.2 2022 October 11 02:33:30 ^ 0.782 < 5.24 x 107° — 
20222840444 142.0 2022 October 11 04:44:25 0.362 < 8.01 x 107? — 
20222840548 145.8 2022 October 11 05:47:44 1.499 < 3.98 x 107° — 
20222840735 152.3 2022 October 11 07:34:38 1.500 < 4.04 x 107? = 
20222840915 158.3 2022 October 11 09:14:53 — 1.254 < 4.59 x 107? — 
20222841047 163.8 2022 October 11 10:47:19 1.453 < 4.26 x 10-9 — 
00015314126 20222840250 135.2 2022 October 11 02:50:06 0.806 < 7.35 x 10-9 — 
00015314127 20222840858 157.3 2022 October 11 08:58:33 0.899 < 6.83 x 107° — 


Table 4. Swift BAT Daily and Total Mosaics of GRB 221009A. Upper limits are 
calculated with an assumed photon index of 1. 
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Start Time Bin End Time Bin Exposure 14-195keV Flux Photon Index 
(ks) (erg cm"? s7!) 

2022 October 08 2022 October 09 0.661 « 2.63 x 107? — 

2022 October 09 2022 October 10 — 6.541 . 5.05 1050 x10 ?  2.30*00» 

2022 October 10 2022 October 11 7.710 4 6.68*21 x 1077° au 

2022 October 11 2022 October 12 4.136 < 6.47 x 107? — 

2022 October 08 2022 October 12 — 19.048 — 1.907077 x 107? — 228*015 


To + 0.046 days 
To + 0.253 days 
To + 0.370 days 
To + 0.570 days 


my 
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Figure 7. XRT WT-mode 1-D PSF profiles as a function of time in the 0.8-5.0 keV band, illustrating the evolution of the echo 
in the early WT data. 


—4 


height of the XRT CCD, but only the central 200 columns are read out, which can make a suitable background region 
difficult to obtain. 

Example time-sliced 1-D profiles from the W'T-mode data are shown in Figure 7. These suggest that within a 
20-pixel source accumulation region, the data are dominated by the source. The profile from the first snapshot of WT 
data is shown in more detail in Figure 8, along with a fit of the profile expected from a point source; the lower panel 
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Figure 8. Upper panel: XRT WT mode 1-D PSF profile from the first snapshot of data at To + 3.4ks (solid line) and the 
modelled point source profile (dashed line). Lower panel: the data minus model residual. 


shows the difference between the data and model. As the GRB was so bright in this snapshot, the central core profile 
is suppressed due to pile-up (resulting in a negative residual, which is not shown for clarity); the central 10 columns 
were subsequently excluded in the spectral extraction to remove the piled-up core. The remaining residuals suggest 
the scattering echoes contribute ~ 10 — 15 96 of the observed count rate in a 5 — 20 pixel radius extraction region 
centred on the source; thus the GRB spectral normalization will be overestimated by a similar amount due to the echo 
contamination in this snapshot. 

To investigate the effect of the echo on the subsequent W'T-mode flux measurements, we took the PC mode image 
from Tp + 1.036 — 1.111 day and summed it in one dimension so as to mimic the WT mode profiles. We did this 
twice; first using all data, then a second time after removing a model of the central source PSF, which was fit to the 
inner-most 20-pixel radius data. The 1-D intensity profiles thus obtained show the source-less data are approximately 
flat, with a variation of « 25 96 between the source outer extraction radius and the wings out to 6 arcminutes. Within 
the source accumulation region, ~ 75 96 of the events come from the source. Thus, the source light curve is not 
significantly altered by the variation in the estimated background level obtained from the wings of the 1-D profile 
where the dust dominates. 
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Table 5. An estimate of the impact of dust on the W'T-mode spectral fit results; see text for details. 
1 (PC-WT)/PC 


Parameter PC result Pseudo-WT result Percentage Difference! 
Nu (10?? cm~?) 1.1-- 0.2 1.43 0.2 27% 
Photon index 1.47 + 0.08 1.57 + 0.09 7% 
Flux (0.3-10keV, 107*° erg cm? s^!) 1.64 + 0.06 1.74 + 0.06 6% 


HER — 


-0.95 -0.65 0.54 5.2 24 


Figure 9. The extracted PC-mode data designed to mimic the WT-mode extraction to investigate the impact of dust. Left: 
the source region, right: the background region. The regions are tilted at the roll angle of the spacecraft. 


To determine the impact of dust on the W'T-mode spectroscopic data, we used the PC mode data in ObsID 
01126853004, comprising 3.1 ks of data collected between 89 and 101 ks after the Fermi trigger.’ First we extracted 
data from a vertical region corresponding to a typical W'T-mode background region, and also from a 20-column radius 
region around the source, excluding the source itself. The spectrum of the dust in the source columns is visibly slightly 
harder than that in the background region; fitting an absorbed power-law in XSPEC showed the difference to be 
minimal, however, with the best-fitting parameters varying by less than 696 between the two spectra and comfortably 
agreeing within their 90%-confidence errors. To better quantify the effect of this possible change in dust spectrum 
between source and background regions, we extracted two spectra from the PC mode data. The first was taken from 
a 20-pixel radius circle centred on the source (i.e. with no dust present) with a background spectrum taken from the 
dust-free region identified above. In the second case we mimicked WT data, extracting data from the full vertical 
height of the CCD for regions corresponding to those used in the WT-mode analysis (Figure 9). We fit these two 
spectra independently in XSPEC using the absorbed power-law model. The results are shown in Table 5 and in 
Figure 10. While all parameters agree to within their errors, the best-fit photon index is 7% higher (i.e., softer) and 
the absorption column 2796 higher in the WT-style data. Since these data were gathered when the source was fainter 
than in the real WT data (i.e. the impact of dust contamination is at its worst in this experiment), we adopt these 
percentages as the maximum inaccuracy expected from our W'T-mode spectral fits. 


7 Although the first snapshot from this ObsID started at 68 ks, its 
ontime was only 7.5 s, meaning the exposure effectively started 
during the second snapshot at 89 ks post trigger. 
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Figure 10. A comparison of a PC-mode spectrum (black), in which dust has been correctly accounted for, and a mimicked 
W'T-mode spectrum (red) of the same data, enabling us to quantify the impact of the dust on the WT mode data. 
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Figure 11. GSC 2-20keV images obtained with GSC_5 at To--2.5ks (left) and GSC.4 at 7)+13.6ks (right). Source and 
background regions are shown by the solid lines and the dashed lines, respectively. 


C. UVOT AFTERGLOW PHOTOMETRY 
The final photometry measured for GRB 2210094 is displayed in Table 8. 


D. MAXI SPECTRAL ANALYSIS AND DUST SCATTERING MITIGATION 


Source spectra at the first (754-2.5 ks) and second (75-4-8.0 ks) scans were extracted from a 3.0? x 4.0? rectangular 
region centered on the source (Figure 11, left), corresponding to FwhiteM of 1.5 deg and 1.5-2 deg of the PSF for 
the scan and its perpendicular (anode) directions, respectively (Sugizaki et al. 2011). Since the source was detected 
near the center of the GSC_4 detector and at the edge of the GSC. 5 detector (f = 2-3 deg, see Figure 2 in Mihara 
et al. 2011), we extracted the background spectra from two 2.4? x 4.0? rectangular regions before and after scanning 


Table 6. MAXI observation logs and spectral fit results. 
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T — To Flux (1078 erg em ? s^!) Dust flux 
Photon index C-stat /d.o.f 
(ks) 0.3-10 keV | unabs 0.3-10 keV — 4-10 keV 0.3-10 keV 
Single power-law fit (Range : 2-20 keV) 
2.459 1.95+0.05 — 6.98 +0.21 ÍBA 3.68 + 0.11 — 353/367 
1.857 6.83 + 0.19 12.32 + 0.34 SA REIS — 356/368 
8.033 RA 1.31298 BOT 0.60 + 0.05 — 180/208 
1.85f 19955 2:21 + 0.15 0.68 0.05 — 187/209 
13.608 1.85f 0.66 + 0.08 1.19+9:15 26 — 93/84 
19.181 1.857 0.45 + 0.07 0 0.25 + 0.04 — 59/62 
24.757 1.857 0.37+9:07 0.06704 G20") as — 40/51 
30.330 1.857 0.12 + 0.06 2 0.07 + 0.03 — 54/37 
Single power-law fit (Range : 4-20 keV) 
2459  1.75+0.09 (01731 IP: 3.46 + 0.13 — 307/322 
1.85/ 6.34 + 0.23 11.4419:42 3.50022 — 309/323 
8.033 0 ee bip 2.597 24 0.58 + 0.06 — 145/166 
1.857 1.05 + 0.10 207, 0.58 + 0.06 — 146/167 
13.608 1.857 0.60 + 0.12 1.08" 534 0:33" 5 a — 42/36 
19.181 1.857 0.36+9:99 0.65+9-16 0.207505 — 43/41 
24.757 1.857 0.30 + 0.08 0.55*5 11 017190 — 29/33 
30.330 1.85f UNE M 0.244918 0.07* 5:04 — 52/32 
Two power-law fit 
2.459 1.81+0.12 6.167983 10.791766 SAT TS liste 351/366 
1.857 | 0.367025 11.48*029 3-510355 0.914937 351/367 
8.033 BA Eh A se PR 059 5 UN Be 180/207 
1.857 0.98 + 0.13 1.77 + 0.23 0.54 + 0.07 ü.4o* 9-20 181/208 
13.608 1.857 ÚT 005 0.29 + 0.07 0.24+9:29 92/83 
19.181 1.857 MT 5 0.49+9:19 0.151008 0.321917 55/61 
24.757 1.857 0.25+9:10 tat 0.14 + 0.05 0.224917 38/50 
30.330 1.857 0.12 + 0.06 0.22 + 0.11 0.07+9:93 0.00 + 0.00 54/36 


NoTE—Flux columns show observed flux at 0.3-10 keV, unabsorbed flux at 0.3-10keV, and observed 
flux at 4-10 keV. All the exposure times are 47s except for 48s at 21:42:30. All the fluxes are in units 


of 10 5 erg em * s^!. The photon index of 1.85 with no error is a fixed value. Note 
are lo. 


the source region, avoiding a shadow region near the center frame of the GSC_4 camera bo 


all errors given 


dy (at 8 = 0) and a high 


background region in the GSC_5. GSC_4 spectra at or after the scan at 17:04 (To4-13.6 ks) were obtained from a 
circular region within a radius of 1.5 arc-deg, and background ones were from an annulus region with the inner and 
outer radii of 2.0 and 4.0 arc-deg overlapped with a 8.0? x 3.4? rectangular region (Figure 11, right). GSC.5 spectra 


at those scans were not used because of low signal-to-noise data. 
We evaluated fluxes of the direct afterglow component in the following two ways: First, w 


e fitted the spectra in the 


4-20 keV band where the contribution from the dust scattering component is small (middle rows in Table 6). This led 
to harder fitted photon indices in the first and second scan spectra, I' = 1.75 + 0.09 and T = Sort. than those in 


the simple 2-20 keV fits, the former being closer to that of the first XRT observation, 1.61 4 


- 0.02. 


Next we directly investigated energy spectra of the dust scattering component using the XRT data. To emulate the 


dust scattering spectrum, we fit the XRT spectra of two bright outer dust-echo rings at To 


+95.4ks and To+146.7 ks 
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Figure 12. GSC spectra at 754-2.5 ks (left panel) and at To+8.0ks (right). The blue solid lines show the direct power-law 
component with T free for the spectrum at To+2.5 ks and fixed at 1.85 at To+8.0 ks. The light blue dashed lines are the power- 
law model with T = 3.94. The red points and lines are background data and models, respectively. The background model is the 
sum of two power-laws. The best-fit data to the model ratio is also shown in each panel. 


with the above power-law model with the fixed column densities. The resultant weighted mean power-law spectral 
index is aust = 3.94 + 0.04. Assuming this value does not change with time, we attempted to fit the GSC spectra at 
2-20 keV with a model with two power-laws, allowing the normalization of both components and the non-dust photon 
index for the first and second scan spectra to be free (Figure 12). We summarize the fitting results in the lower rows 
in Table 6. Interestingly, the obtained parameters for the first and second scan spectra are almost perfectly consistent 
with those in the single power-law fits at 4-20 keV even if the photon index is fixed or not, though the uncertainty of 
the dust flux (the normalization of the power-law with T = 3.94) is large. We also note that the first scan spectrum 
tends to be harder than the second scan one in all the cases. 

The differential scattering cross-section has energy E and scattering angle 0, dependence of approximately 
exp(-aF?0?) for 6; < 1 where a is a constant of proportional to the square of the size of grain (e.g., Mauche & 
Gorenstein 1986). Thus, an actual dust scattered component is considered to have a harder spectrum than assumed 
here especially in the first and second scan observations. If the time and spectral evolution of the dust echo ring 
is understood, the fitting parameters can be more constrained. This is, however, beyond the scope of this paper. 
Finally, we note that if we assume I = 3.5 for the scattered power-law component we got a steeper photon index 
(T = 1.75 + 0.14) and a lower 0.3-10 keV absorbed flux of 5.7619:32 x 1078 erg em ? s~! for the direct power-law 
component in the first scan spectrum. 


E. DUST ECHO MODELING 


In order to examine the properties of the dust scattered echo, images were created from the XRT PC mode event 
data over the 0.8 — 5.0 keV band, where the dust reprocessing cross-section peaks. Given the scattered emission evolves 
radially with time, the images were initially extracted over per snapshot intervals (with typical exposures ranging from 
0.3 to 1.6ks), until day 22 post trigger when per ObsID images were created (with exposures from 3.3 to 5.2 ks). 
Vignetting-corrected exposure maps were also generated, spanning the observing times for each image. 

Radial profiles were then generated as follows. The GRB position was obtained on a per image basis by fitting the 
expected XRT PSF model (including a piled-up profile modification when required, as described in Evans et al. 2020) 
to the imaging data out to 47.1" where the central source dominates. This position was then used as the location 
about which radial profiles were generated, initially in 2.357" linear bins from 0.79 to 13.75 ' radius, following the 
removal of 12 faint point sources (with count rates less than 3.5 x 10^? count s~!, determined from 67 ks of late time 
data ranging from To 4- 43 to 64 days after the trigger). The central source profile was retained and combined with 
a nominal background level determined from late time data to provide the non-halo background estimate. Example 
radial profiles are shown in Figure 13 which shows the echo emission expanding and decreasing in intensity with time. 
Later profiles were stacked to increase signal to noise in the plot, with radial bins scaled so that radial bins correspond 
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Figure 13. XRT PC mode radial intensity profiles in the 0.8-5.0 keV band (blue solid lines); later profiles show stacked 
observations in indicated time windows, with radii scaled to the beginning of the observing window using eq. (1) so that radial 
intensity features at the same distance appear at the same radius. Overplotted lines show the PSF and background model 


(green dashed) and the dust model (orange dashed) and the best fit model (solid red). 


to the first profile in the set of stacked observations according to eq. (1) in order to align radial intensity features of 


different observations. 


The plots show clear excess over the PSF from the afterglow emission by the point source, with each peak (i.e., each 


ring in the image) corresponding to a separate dust cloud along the line of sight. 


Approximating the light curve of the prompt emission as a delta function, the radial intensity of the echo at time 
t = terg + At plotted in Figure 13 can be written as 


Techo = 


2c- n(DDust) do Lg 
02 . F . dQ . € 
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where n(D) is the dust volume density as a function of distance D(0, At) along the line of sight for a given scattering 
angle 0 and time since the burst At and Tphot is the optical depth to photo-electric absorption along the line of sight. 
The observed radial intensity distribution is then given by the convolution of eq. (E1) with the telescope PSF. 

The flux of a given cloud of column density N at distance D (used to perform the dist model fits) is given by (Heinz 

et al. 2016) 
2ncN do 
lone ^3 

Because the time delay from the prompt emission and the radius 0 of each ring are known, eq. (E1) can be used to 
measure the product of burst fluence, density, and scattering cross section, and, with knowledge of the fluence and a 
model for the scattering cross section, we can, in principle, solve for the dust density distribution nanst as a function 
of line-of-sight (LOS) distance and the azimuthal angle along the ring for each observation. 

In practice, estimating the fluence at soft X-ray energies (where dust scattering is most efficient and where the echo 
is observed) requires extrapolation from the sensitivity band of the instruments that detected the prompt emission 
(namely, GBM; Lesage et al. in prep.), corrected for intrinsic absorption. This fluence value, derived using GBM data 
(after pile-up correction procedures), has an uncertainty of — 5096 that propagates to our results. 

In addition, models for the dust scattering cross section vary from cloud to cloud and are themselves not well 
constrained. Thus, a conservative approach will treat both do/dQ and npust as functions to be fitted simultaneously, 
with the understanding that constraints derived for each will carry some degeneracy. 

'That said, the relative column density distribution is well constrained by this process, provided that dust chemistry 
and grain size distribution do not vary drastically from cloud to cloud. 

At typical soft X-ray energies around 2keV, the scattering cross section is roughly constant at angles smaller than 
about 100" while falling off rapidly at larger scattering angles (roughly as do/d€) œ 935). This results in a rapid 
decline in the intensity and flux of the echo as a function of time as the ring size expands, as can be seen in Figure 13. 
It further implies that, beyond a ring radius of about 100", smaller rings, which are farther away, will have a larger 
intensity per column density compared to larger rings, which are closer to the observer. 

On the other hand, the Galactic latitude of GRB 221009A of b = 4.3? suggests that the Galactic column density 
distribution should be dominated by nearby dust within the plane, within a few kpc. For example, dust at the 
intersection of the LOS with the Solar circle on the far side of the Galaxy would lie 780pc (almost 8 dust scale heights, 
e.g., Li et al. 2018) above the Galactic plane, where we would not expect to see substantial amounts of dust (however, 
any dust present at those distances will benefit from the intensity enhancement discussed in the previous paragraph). 

'To constrain the dust column density and cross section models, we decomposed the dust distribution along the 
line of sight into 100 logarithmically spaced bins and determined the least squares fit for the column density for each 
bin. To this end, we constructed 1-D XRT point spread functions to approximate the radial profile for each dust 
ring. We then defined the fitting function as given by eq. (E1) for each column density component. We added the 
point source PSFs for each profile to the fitting function. We then fitted all 100 components simultaneously to all the 
background-subtracted and radial intensity profiles between day 1 and 42 post-burst. The resulting inferred density 
distribution is plotted in Figure 4. 


(E2) 


F. BROADBAND SED CONSTRUCTION 


The BAT survey spectra are produced using the BatAnalysis python package (Parsotan et al, in prep.). Using 
the HEASoft batsurvey script, the BatAnalysis package calculates the count rate of the GRB in each of the 8 
energy bins used in the BAT survey (14-20, 20-24, 24-35, 35-50, 50—75, 75-100, 100-150, and 150-195 keV) and the 
errors associated with each energy bin. The package also generates the detector response matrix using the HEASoft 
batdrmgen script for the BAT survey observation of interest. 

For the time intervals of the two SEDs, we created time-sliced X-ray spectra being careful to limit dust echo 
contamination (see 83.2). The source spectral files were grouped to > 20 counts per energy bin. The spectral files were 
normalized to correspond to the 0.3-10 keV flux of the afterglow at To 4- 4.2 ks or To +43 ks. The flux used to normalize 
a given spectrum is determined by fitting a power-law to the data within the SED time range and the best-fit decay 
index is used to compute the flux at the mid-point of the SED, in the same way as was done for the UVOT data. 

To build the optical spectral files we follow the methodology provided in Schady et al. (2010). This essentially sets 
the count rates in the spectral files to that determined at an instantaneous epoch by extrapolating the UVOT light 
curves. In order to obtain the count rates at the instantaneous epoch, we first combined the data from the different 
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Table 7. Detectable redshift limit for each instrumental setup. “Grid ID" indicates the location at the BAT image plane and 
corresponds to different partial coding fractions and thus different detector sensitivities. ^ndet" refers to the number of enabled 
detectors. 


Grid ID ndet Zim 


17 29000 11 
17 22000 11 
17 18000 11 
17 15000 11 
16 29000 11 
16 22000 11 
16 18000 11 
16 15000 11 
15 29000 11 
15 22000 10 
15 18000 9 
15 15000 9 
14 29000 5 
14 22000 4 
14 18000 4 
14 15000 3 


filters into a single-flow light curve filter. The light curves of the different filters were normalized to that in the v 
filter. The normalization was determined by fitting a power-law to each of the light curves in a given time range 
simultaneously. The power-law indices were constrained to be the same for all the filters and the normalizations were 
allowed to vary between the filters. The ratios of the power-law normalizations between each filter and the v filter 
were then used to shift the individual light curves to the same normalization as the v filter, thus resulting in a single 
filter light curve. To construct the two SEDs at Tp + 4.2 ks and Ty + 43 ks, we first determine the temporal slope from 
the single filter light curve within the corresponding time interval. By fixing the power-law index at this value, we 
then fit a power-law to the individual filter light curves. We use the derived normalizations to compute the count rate 
and count rate error at the required time, which was then applied to the relevant spectral file. 


G. BAT TRIGGER SIMULATION 


To estimate the intrinsic rate of such events, we utilize the BAT trigger simulator (Lien et al. 2014) to calculate 
the detectability of the prompt emission at different redshifts under different representative instrumental setups. 
Specifically, our setup uses (1) one standard average background level, (2) four different locations on the BAT image 
plane (i.e., Grid ID of [14, 15, 16, 17] or equivalent boresight angle of [56,45,27,0] deg), which represent different 
detector sensitivity and cover different locations within the BAT field of view, and (3) four different numbers of 
enabled detectors (29000, 22000, 18000, and 15000; these numbers represent the change of average number of enabled 
detectors from 2005 to 2022). We ran simulations with a combination of each of these setups with a sample of redshifts 
from z = 0.1 to z = 12 with increments of Az = 1.0 (i.e., z = 0.1,1.0,2.0,...,12.0). The simulation results are 
summarized in Table 7. 

For each setup, we estimate the expected number of detections based on the highest detectable redshift zm and an 
assumed intrinsic rate with the following equation: 


Näet isetup = Rana (z < Zlim) Ísurvey Jio feria Índet- (G3) 


Rerp(Z < zin) is the all-sky intrinsic rate up to the redshift limit for this burst to be detected by BAT with a specific 
setup of Grid ID and number of enabled detectors, and is calculated by integrating the comoving rate by taking into 
account of the volume of the universe through the following equation (see, e.g., Lien et al. 2011, for details of the 
derivation): 

2 


Zlimit re dr 
R im) = Í R comov(Z = comov “! comov " 
GRB(Z < Zim) T f GRB, (z=0) a dz (G4) 
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where the comoving distance reomov(z) = ( 5 ) p l/H(z)dz. fsurvey is the fraction of time that BAT is capable 
of triggering, and we adopt fsurvey = 0.8 based on the study in Lien et al. (2016). fr, = (2.1sr)/(4msr) is the 
fraction of sky covered by the entire BAT field of view down to a partial coding fraction of = 0.1, and feria = 
[0.460, 0.346, 0.109, 0.085] is the fraction of BAT field of view for Grid ID 14, 15, 16, and 17 8. In other words, 
Ítov X feria is the fraction of bursts in the entire sky that would have the partial coding fraction of the specific Grid 
ID. fndet is the fraction of GRBs in the BAT field of view that would occur with this number of enabled detectors. 
For simplicity, we assume an equal number of GRBs occur with these four numbers of enabled detectors. That is, 
fadet = 0.25. 

The total number of detections can then be calculated by adding up Naet,isesup for all 16 combinations of instrumental 
setups listed in Table 7. That is, 


16 

Naet,tot = 5 Naeti ccs (G5) 
isetup =1 

We assume a flat intrinsic comoving rate of RerB,comov(z = 0) and adjust the value until the detection rate matches 
with that of a GRB 221009A-like event. We set the upper limit of the BAT detection rate of such events to be 1 per 
18 years of the Swift mission lifetime, because the prompt emission of the burst actually occurred outside of the BAT 
field of view. This gives us a corresponding upper limit of RarB.comov(z = 0) € 6.1 x 1074 Gpe~3 yr7!. Integrating 
this flat comoving rate from z = 0 to z = 12, we obtain an upper-limit on the all-sky rate of GRB 221009A-like events 
to be 0.5 yr Í. Comparing to the all-sky intrinsic long-GRB rate of ~ 4571 yr~! in Lien et al. (2014), the fraction of 

GRB 221009A-like events is roughly 0.5/4571 < 1.0 x 1074. 


Table 8. Swift UVOT observations. 


Tmia — To Half Exposure Magnitude Flux Filter S/N 
(ks) (s) (erg em ? s^! ÅT!) 
3.4 25 16.743*0044 9.64-E0.44 x 1071 white 21.956 
3.5 30 16.762*0049 9.48 +0.40 x 10716 white 23.877 
3.5 20 16.726+9:957 9.80 +0.50 x 10716 white 19.698 
3.8 10 16.714085 9.94 +0.85 x 1071 white 11.664 
4.0 10 16.95810'094 7.912: 0.72 x 1016 white 11.024 
4.1 75 16.972*0075 7.812: 0.52 x 1019 white 15.078 
4.4 10 16.9076*0095 7.78 +0.71 x 10716 white 10.976 
21.9 246 18.803*006? 1.45 +0.08 x 1071" white 17.686 
44.9 76 19.862*0199 5.45 +0.77 x 1071" white — 7.065 
61.2 408 20.120053 4.38 +0.33 x 10!" white 13.261 
120.0 14394 21.152*007? 1.66 £0.12 x 107!" white 14.176 
152.3 17859 21.655*0099  1.0523:0.09 x 107! white 11.522 
198.5 23982 21.599*072. 1.10+0.19x 1071” white — 5.685 
258.6 31606 21.987019,  7.712- 1.24 x 10 75. white — 6.218 
306.9 6602 22.574*0322 4491.16 x 101? white 3.872 
436.1 49190 22.703*0 259 3.99 1.44 x 10775 white 2.773 
585.3 71725 22.6881036 4.04 1.59 x 1071? white 2.541 
754.3 91946 22.504*0272 4.41 +1.81 x 1077? white 2.440 


Table 8 continued 


8 https://asd.gsfc.nasa.gov/Craig.Markwardt/bat-cal/ 
solid-angle/ 


Table 8 (continued) 


Tmia — To Half Exposure Magnitude Flux Filter S/N 
(ks) (s) (erg em ?s TÁ 
971.9 120245 > 22.625 «4.28 x 10718 white | — 
1251.6 154713 22.882 «3.38 x 10715 white | — 
1701.2 192141 > 22.836 <3.52 x 10775 white | — 
2087.8 137716 > 23.166 «2.60 x 10715 white | — 
3940.8 432429 > 23.602 «1.74 x 10718 white | — 
3.4 5 15.463*01€92 2.44 + 0.34 x 107" v 7.210 
3.9 10 15.53819:123 2.28 +0.25 x 107" v 9.220 
4.2 114 15.720'0 peg” 1.91 +0.16 x 107" v 12.072 
4.4 10 15.893*0129 1.64 0.21 x 107" v 7.995 
38.2 412 18.152 19 0:5, 2.05 +0.14 x 10716 v 14.642 
55.3 412 18.755191! 1.18 0.11 x 19775 v 10.375 
92.8 8959 19.427012 6.34 0.68 x 1071” v 9.342 
436.6 49188 > 20.443 «2.49 x 10777 v — 
573.7 59918 21.082*099* 1.88 + 0.66 x 10-1” v 2.109 
822.1 97374 > 20.64 «2.08 x 1071" v — 
1065.7 122800 > 21.411 «1.02 x 10777 v — 
1372.4 160437 > 21.617 <8.44 x 10718 v — 
1833.6 220177 > 22.032 «5.76 x 10715 v = 
2151.2 74786 > 21.125 «1.33 x 1077 v — 
3940.0 432472 > 21.862 «6.74 x 10715 v — 
3.8 10 17.062*0122 9.70 Æ+ 1.20 x 10716 b 8.058 
4.0 10 17.22740- 138 8.34 1.11 x 10777 b 7.495 
4.4 92 17.5511915% 6.18 +0.77 x 10716 b 8.073 
44.4 454 20.268 133 5.06 0.63 x 10 7 b 7.999 
60.3 453 20.666+9:237 3.51 + 0.69 x 197" b 5.100 
438.9 52551 > 21.71 «1.34 x 10717 b — 
584.9 71610 > 21.597 «1.49 x 10777 b — 
753.8 92017 22.4441031% 6.82 + 2.56 x 10718 b 2.662 
971.7 120321 > 22.02 «1.01 x 10717 b — 
1251.2 154869 > 22.124 «9.16 x 10715 b — 
1700.9 192324 22.229 «8.32 x 10718 b — 
2087.2 137713 > 22.48 «6.60 x 10715 b = 
3849.4 438307 > 22.969 «4.21 x 10715 b — 
3.6 38 1559210115. 3.24 + 0.36 X 10^ 79 u 8.927 
3.7 50 17.67319113 3.01 Æ 0.30 x 10716 u 9.933 
3.7 37 17.629*0131 3.13 Æ 0.36 x 107° u 8.773 
3.9 10 17.882+9:309 2.48 + 0.61 x 10716 u 4.041 
4.4 97 17.765*0198 2.77 +0.46 x 10-18 u 6.043 
32.4 415 20.258 10 727 2,78 + 0.48 x 10 17 5.769 
49.7 414 20.729*0313  1.80-E0.45 x 1071" u 3.992 


Table 8 continued 
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Table 8 (continued) 


Tmia — To Half Exposure Magnitude Flux Filter S/N 

(ks) (s) (erg cm~? s7! À-1) 

95.8 414 21305609 1.27 + 0.42 x 1077 u 3.010 
438.8 52496 > 20.903 «1.54 x 10717 u = 
584.7 71385 > 20.647 «1.94 x 10717 u — 
753.5 91837 > 21.678 <7.52 x 10718 u — 
971.4 120396 > 21.285 «1.08 x 10717 u = 
1251.1 154783 > 21.32 «1.05 x 1071” u = 
1700.8 192285 > 21.335 «1.03 x 10717 u — 
2086.9 137473 20.098 «3.23 x 1071" u - 
3848.7 438428 22.798 «2.68 x 10718 u == 

3.9 10 18.2381967} 2.02 +0.93 x 10779 — uvw1 2.168 

4.4 97 > 19.108 «9.06 x 1071” uvul — 

31.6 450 > 20.175 «3.39 x 1071” uvul = 

48.8 450 > 20.003 «3.97 x 1071" uvul = 

67.5 450 > 20.5 «2.51 x 10717 uvul = 

92.4 2910 > 21.471 «1.03 x 1071” uvul = 
685.4 269 > 20.833 «1.85 x 1071” uvw1 | — 
1458.7 22863 > 22.436 «4.23 x 10718 uw — 

3.9 10 > 18.228 «2.37 x 10718 uv m2 | — 

4.4 97 > 18.771 «1.44 x 10718 uv m2 | — 

87.1 2398 > 21.314 «1.38 x 1071" uvm2 | — 

3.9 10 18.852 «1.54 x 10716 uvw2 — 

4.1 114 19.356 «9.68 x 1071” uw  — 

4.4 10 18.966 «1.39 x 107 !8 uw  — 

25.6 83 > 19.614 «1.64 x 10717 uw  — 

37.3 450 > 20.917 «2.30 x 1071” uw  — 

54.4 450 > 20.68 «2.86 x 1071” uvw2  — 

73.3 253 > 21.153 «1.85 x 1071” uw  — 
100.5 450 20.525 «3.30 x 1071” uvw2 — 


NoTE—Columns 1 and 2 give the mid-time of the exposure in kiloseconds since the GBM 
trigger and the length of the exposure divided by 2. Magnitudes are given in the Vega system. 
Uncertainties are given at the 1o level. Observations with S/N greater than or equal to 2 were 
considered detections. For non-detections with S/N less than 2, 30 upper limits are given. 
'The values in this table have not been corrected for Galactic extinction. 


REFERENCES 


Arnaud, K. A. 1996, in Astronomical Society of the Pacific Astropy Collaboration, Robitaille, T. P., Tollerud, E. J., 
et al. 2013, A&A, 558, A33, 
doi: 10.1051/0004-6361/201322068 
Barthelmy, S. D., Barbier, L. M., Cummings, J. R., et al. 
17 2005, SSRv, 120, 143, doi: 10.1007/s11214-005-5096-3 


Conference Series, Vol. 101, Astronomical Data Analysis 


Software and Systems V, ed. G. H. Jacoby & J. Barnes, 


Beardmore, A. P., Willingale, R., Kuulkers, E., et al. 2016, 
MNRAS, 462, 1847, doi: 10.1093/mnras/stw1753 

Bloom, J. S., Frail, D. A., & Sari, R. 2001, AJ, 121, 2879, 
doi: 10.1086/321093 

Bloom, J. S., Perley, D. A., Li, W., et al. 2009, ApJ, 691, 
723, doi: 10.1088/0004-637X/691/1/7723 

Breeveld, A. A., Landsman, W., Holland, S. T., et al. 2011, 
in American Institute of Physics Conference Series, Vol. 
1358, American Institute of Physics Conference Series, 
ed. J. E. McEnery, J. L. Racusin, & N. Gehrels, 373-376, 
doi: 10.1063/1.3621807 

Brethauer, D., Grefenstette, B., Racusin, J., et al. 2022, 
GRB Coordinates Network, 32695, 1 

Burrows, D. N., Hill, J. E., Nousek, J. A., et al. 2005, 
SSRv, 120, 165, doi: 10.1007/s11214-005-5097-2 

Cash, W. 1979, ApJ, 228, 939, doi: 10.1086/156922 

Castro-Tirado, A. J., Sanchez-Ramirez, R., Hu, Y. D., 
et al. 2022, GRB Coordinates Network, 32648, 1 

Cenko, S. B., Frail, D. A., Harrison, F. A., et al. 2011, ApJ, 
732, 29, doi: 10.1088/0004-637X /132/1/29 

D'Avanzo, P., Ferro, M., Brivio, R., et al. 2022, GRB 
Coordinates Network, 32755, 1 

de Pasquale, M., Evans, P., Oates, S., et al. 2009, MNRAS, 
392, 153, doi: 10.1111/j.1365-2966.2008.13990.x 

de Ugarte Postigo, A., Izzo, L., Pugliese, G., et al. 2022, 
GRB Coordinates Network, 32648, 1 

Dichiara, S., Gropp, J. D., Kennea, J. A., et al. 2022, GRB 
Coordinates Network, 32632, 1 

Evans, P. A., Beardmore, A. P., Page, K. L., et al. 2007, 
A&A, 469, 379, doi: 10.1051/0004-6361:20077530 

—. 2009, MNRAS, 397, 1177, 
doi: 10.1111/j.1365-2966.2009.14913.x 

Evans, P. A., Page, K. L., Osborne, J. P., et al. 2020, ApJS, 
247, 54, doi: 10.3847/1538-4365 /ab7db9 

Foreman-Mackey, D., Farr, W., Sinha, M., et al. 2019, The 
Journal of Open Source Software, 4, 1864, 
doi: 10.21105/joss.01864 

Frail, D. A., Kulkarni, S. R., Sari, R., et al. 2001, ApJL, 
562, L55, doi: 10.1086/338119 

Galama, T. J., Vreeswijk, P. M., van Paradijs, J., et al. 
1998, Nature, 395, 670, doi: 10.1038/27150 

Gehrels, N., Chincarini, G., Giommi, P., et al. 2004, ApJ, 
611, 1005, doi: 10.1086/422091 

Gendreau, K. C., Arzoumanian, Z., & Okajima, T. 2012, in 
Society of Photo-Optical Instrumentation Engineers 
(SPIE) Conference Series, Vol. 8443, Space Telescopes 
and Instrumentation 2012: Ultraviolet to Gamma Ray, 
ed. T. Takahashi, S. S. Murray, & J.-W. A. den Herder, 
844313, doi: 10.1117/12.926396 


29 


Ghisellini, G., Ghirlanda, G., Nava, L., & Firmani, C. 2007, 
ApJL, 658, L75, doi: 10.1086/515570 

Guha, A., & Nicholson, P. 2022, GRB Coordinates 
Network, 32745, 1 

Heinz, S., Corrales, L., Smith, R., et al. 2016, ApJ, 825, 15, 
doi: 10.3847 /0004-637X/825/1/15 

Heinz, S., Burton, M., Braiding, C., et al. 2015, ApJ, 806, 
265, doi: 10.1088/0004-637X/806/2/265 

Huang, Y., Hu, S., Chen, S., et al. 2022, GRB Coordinates 
Network, 32677, 1 

Iwakiri, W., Jaisawal, G. K., Younes, G., et al. 2022, GRB 
Coordinates Network, 32694, 1 

Izzo, L., Saccardi, A., Fynbo, J. P. U., et al. 2022, GRB 
Coordinates Network, 32765, 1 

Kann, D. A., & Agui Fernandez, J. F. 2022, GRB 
Coordinates Network, 32762, 1 

Kennea, J. A., Williams, M., & Swift Team. 2022, GRB 
Coordinates Network, 32635, 1 

Kobayashi, K., Negoro, H., Nakajima, M., et al. 2022, GRB 
Coordinates Network, 32756, 1 

Kuin, N. P. M., Dichiara, S., & Swift/UVOT Team. 2022, 
GRB Coordinates Network, 32656, 1 

Kulkarni, S. R., Frail, D. A., Wieringa, M. H., et al. 1998, 
Nature, 395, 663, doi: 10.1038/27139 

Laha, S., Wadiasingh, Z., Parsotan, T., et al. 2022, ApJ, 
929, 173, doi: 10.3847/1538-4357 /acbf3c 

Li, L., Shen, S., Hou, J., et al. 2018, ApJ, 858, 75, 
doi: 10.3847/1538-4357 /aabaef 

Lien, A., Chakraborty, N., Fields, B. D., & Kemball, A. 
2011, ApJ, 740, 23, doi: 10.1088/0004-637X/740/1/23 

Lien, A., Sakamoto, T., Gehrels, N., et al. 2014, ApJ, 783, 
24, doi: 10.1088/0004-637X/783/1/24 

Lien, A., Sakamoto, T., Barthelmy, S. D., et al. 2016, ApJ, 
829, 7, doi: 10.3847/0004-637X/829/1/7T 

Liu, Y.-Q., Modjaz, M., Bianco, F. B., & Graur, O. 2016, 
ApJ, 827, 90, doi: 10.3847/0004-637X/827/2/90 

Margutti, R., Milisavljevic, D., Soderberg, A. M., et al. 

2014, ApJ, 797, 107, doi: 10.1088/0004-637X /797/2/107 

Maselli, A., Melandri, A., Nava, L., et al. 2014, Science, 

343, 48, doi: 10.1126/science.1242279 

Mathis, J. S., Rumpl, W., & Nordsieck, K. H. 1977, ApJ, 

217, 425, doi: 10.1086/155591 

Matsuoka, M., Kawasaki, K., Ueno, S., et al. 2009, PASJ, 

61, 999, doi: 10.1093/pasj/61.5.999 

Mauche, C. W., & Gorenstein, P. 1986, Yapj, 302, 371, 

doi: 10.1086/163996 

Meegan, C., Lichti, G., Bhat, P. N., et al. 2009, ApJ, 702, 

791, doi: 10.1088/0004-637X/702/1/791 

Mihara, T., Nakajima, M., Sugizaki, M., et al. 2011, PASJ, 
63, $623, doi: 10.1093/pasj/63.sp3.5623 


30 


Nasa High Energy Astrophysics Science Archive Research 
Center (Heasarc). 2014, HEAsoft: Unified Release of 
FTOOLS and XANADU, Astrophysics Source Code 
Library, record ascl:1408.004. http://ascl.net/1408.004 

Negoro, H., Kohama, M., Serino, M., et al. 2016, PASJ, 68, 
$1, doi: 10.1093/pasj/psw016 

Negoro, H., Nakajima, M., Kobayashi, K., et al. 2022, The 
Astronomer's Telegram, 15651, 1 

Negro, M., Di Lalla, N., Omodei, N., et al. 2023, arXiv 
e-prints, arXiv:2301.01798, 
doi: 10.48550/arXiv.2301.01798 

Nousek, J. A., Kouveliotou, C., Grupe, D., et al. 2006, ApJ, 
642, 389, doi: 10.1086/500724 

Perley, D. A., Cenko, S. B., Corsi, A., et al. 2014, ApJ, 781, 
37, doi: 10.1088/0004-637X/781/1/37 

Piran, T. 2004, Reviews of Modern Physics, 76, 1143, 
doi: 10.1103/RevModPhys.76.1143 

Planck Collaboration, Aghanim, N., Akrami, Y., et al. 
2020, A&A, 641, A6, doi: 10.1051/0004-6361/201833910 

Poole, T. S., Breeveld, A. A., Page, M. J., et al. 2008, 
MNRAS, 383, 627, doi: 10.1111/j.1365-2966.2007.12563.x 

Racusin, J. L., Liang, E. W., Burrows, D. N., et al. 2009, 
ApJ, 698, 43, doi: 10.1088/0004-637X/698/1/43 

Remillard, R. A., Loewenstein, M., Steiner, J. F., et al. 
2022, AJ, 163, 130, doi: 10.3847/1538-3881/ac4ae6 

Roming, P. W. A., Kennedy, T. E., Mason, K. O., et al. 
2005, SSRv, 120, 95, doi: 10.1007 /s11214-005-5095-4 

Sari, R., Piran, T., & Halpern, J. P. 1999, ApJL, 519, L17, 
doi: 10.1086/312109 


Sari, R., Piran, T., & Narayan, R. 1998, ApJL, 497, L17, 
doi: 10.1086/311269 

Schady, P., Mason, K. O., Page, M. J., et al. 2007, 
MNRAS, 377, 273, doi: 10.1111/j.1365-2966.2007.11592.x 

Schady, P., Page, M. J., Oates, S. R., et al. 2010, MNRAS, 
401, 2773, doi: 10.1111/j.1365-2966.2009.15861.x 

Schlafly, E. F., & Finkbeiner, D. P. 2011, ApJ, 737, 103, 
doi: 10.1088/0004-637X /737/2/103 

Schlegel, D. J., Finkbeiner, D. P., & Davis, M. 1998, ApJ, 
500, 525, doi: 10.1086/305772 

Schnoor, P. W., Nicholson, P., & Welch, D. L. 2022, GRB 
Coordinates Network, 32744, 1 

Sugizaki, M., Mihara, T., Serino, M., et al. 2011, PASJ, 63, 
$635, doi: 10.1093/pasj/63.sp3.S635 

Tiengo, A., Pintore, F., Mereghetti, S., & Salvaterra, R. 
2022, The Astronomer's Telegram, 15661, 1 

Vasilopoulos, G., & Petropoulou, M. 2016, MNRAS, 455, 
4426, doi: 10.1093/mnras/stv2605 

Veres, P., Burns, E., Bissaldi, E., et al. 2022, GRB 
Coordinates Network, 32636, 1 

Willingale, R., Starling, R. L. C., Beardmore, A. P., Tanvir, 
N. R., & O’Brien, P. T. 2013, MNRAS, 431, 394, 
doi: 10.1093/mnras/stt175 

Yamazaki, R., Sato, Y., Sakamoto, T., & Serino, M. 2020, 
MNRAS, 494, 5259, doi: 10.1093/mnras/staa1095 

Zhang, B., Fan, Y. Z., Dyks, J., et al. 2006, ApJ, 642, 354, 
doi: 10.1086/500723 


